{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import sys\n",
    "sys.path.append('../GTF3')\n",
    "from Utilities import *\n",
    "from admm import admm\n",
    "\n",
    "from IPython.core.interactiveshell import InteractiveShell\n",
    "InteractiveShell.ast_node_interactivity = \"all\"\n",
    "from functools import *\n",
    "\n",
    "import networkx as nx\n",
    "import pygsp as pg\n",
    "import numpy as np\n",
    "import random\n",
    "from time import time\n",
    "\n",
    "from scipy.sparse import csr_matrix\n",
    "from scipy.linalg import inv\n",
    "\n",
    "from sklearn.model_selection import ParameterGrid\n",
    "from hyperopt import hp\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline\n",
    "from IPython.display import display\n",
    "from math import log\n",
    "#from tqdm import tqdm_notebook, tnrange\n",
    "import pandas as pd\n",
    "import datetime as dt\n",
    "import matplotlib as mpl\n",
    "import matplotlib.pyplot as plt\n",
    "from matplotlib import rc\n",
    "\n",
    "rc('text', usetex=True)\n",
    "mpl.rcParams['text.usetex'] = True\n",
    "mpl.rcParams['text.latex.preamble'] = [r'\\usepackage{amsmath}'] #for \\text command\n",
    "\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    input_snr (dB)  sigma_sq penalty_f  gamma_mult  rho_mult  penalty_param  \\\n",
      "0              -10  112.5000        L1      0.2065    0.1102            0.0   \n",
      "1              -10  112.5000        L1      0.0960    0.1514            0.0   \n",
      "2              -10  112.5000   SCAD-L1      4.5160    1.1985            3.7   \n",
      "3              -10  112.5000   SCAD-L1      7.6041   20.0521            3.7   \n",
      "4              -10  112.5000    MCP-L1     25.0228    0.5140            1.4   \n",
      "5              -10  112.5000    MCP-L1     10.3409   16.7093            1.4   \n",
      "6               -5   35.5756        L1      0.3668    0.1075            0.0   \n",
      "7               -5   35.5756        L1      0.1627    0.3428            0.0   \n",
      "8               -5   35.5756   SCAD-L1      0.6803    0.2695            3.7   \n",
      "9               -5   35.5756   SCAD-L1     15.9216    7.1669            3.7   \n",
      "10              -5   35.5756    MCP-L1      1.5306    0.1525            1.4   \n",
      "11              -5   35.5756    MCP-L1     34.1253    3.5041            1.4   \n",
      "12               0   11.2500        L1      0.6248    5.8824            0.0   \n",
      "13               0   11.2500        L1      0.2881    0.3259            0.0   \n",
      "14               0   11.2500   SCAD-L1      4.5813    0.2867            3.7   \n",
      "15               0   11.2500   SCAD-L1     24.3336    2.0544            3.7   \n",
      "16               0   11.2500    MCP-L1      7.0082    0.1329            1.4   \n",
      "17               0   11.2500    MCP-L1      1.1229    0.2917            1.4   \n",
      "18               5    3.5576        L1      1.2948    4.3710            0.0   \n",
      "19               5    3.5576        L1      0.5383    3.0922            0.0   \n",
      "20               5    3.5576   SCAD-L1      2.7433    0.2019            3.7   \n",
      "21               5    3.5576   SCAD-L1      1.4941    0.7834            3.7   \n",
      "22               5    3.5576    MCP-L1     24.5706    0.8291            1.4   \n",
      "23               5    3.5576    MCP-L1      8.2185    0.4619            1.4   \n",
      "24              10    1.1250        L1      2.1129    4.5473            0.0   \n",
      "25              10    1.1250        L1      0.9405    2.4232            0.0   \n",
      "26              10    1.1250   SCAD-L1     53.9420    2.1853            3.7   \n",
      "27              10    1.1250   SCAD-L1      8.3250    2.0922            3.7   \n",
      "28              10    1.1250    MCP-L1     54.4157    1.1317            1.4   \n",
      "29              10    1.1250    MCP-L1     46.8948    0.4756            1.4   \n",
      "30              15    0.3558        L1      3.8249    6.1959            0.0   \n",
      "31              15    0.3558        L1      1.6889   10.5369            0.0   \n",
      "32              15    0.3558   SCAD-L1     40.0820    5.1895            3.7   \n",
      "33              15    0.3558   SCAD-L1     49.2767    3.1571            3.7   \n",
      "34              15    0.3558    MCP-L1     53.8908    2.0692            1.4   \n",
      "35              15    0.3558    MCP-L1     47.6396    1.2511            1.4   \n",
      "36              20    0.1125        L1      6.9798   14.6076            0.0   \n",
      "37              20    0.1125        L1      2.9636   11.8453            0.0   \n",
      "38              20    0.1125   SCAD-L1     54.2747    8.5047            3.7   \n",
      "39              20    0.1125   SCAD-L1     52.4211    3.7511            3.7   \n",
      "40              20    0.1125    MCP-L1     47.4904    3.7670            1.4   \n",
      "41              20    0.1125    MCP-L1     50.6210    1.4927            1.4   \n",
      "42              25    0.0356        L1     12.1498   14.6964            0.0   \n",
      "43              25    0.0356        L1      5.4368   19.6091            0.0   \n",
      "44              25    0.0356   SCAD-L1     53.8617   15.8310            3.7   \n",
      "45              25    0.0356   SCAD-L1     43.9441   11.4919            3.7   \n",
      "46              25    0.0356    MCP-L1     53.4987    4.7119            1.4   \n",
      "47              25    0.0356    MCP-L1     38.7646    3.3045            1.4   \n",
      "48              30    0.0112        L1     20.1298   19.9748            0.0   \n",
      "49              30    0.0112        L1      9.3838   16.4264            0.0   \n",
      "50              30    0.0112   SCAD-L1     54.2127   14.3469            3.7   \n",
      "51              30    0.0112   SCAD-L1     45.2690   20.0808            3.7   \n",
      "52              30    0.0112    MCP-L1     43.2715   10.9308            1.4   \n",
      "53              30    0.0112    MCP-L1     54.4044    7.7025            1.4   \n",
      "\n",
      "    avg output_snr (dB)        0        1        2        3        4        5  \\\n",
      "0                8.6422   8.3054   8.9744   8.3184   8.8736   7.5636  11.4627   \n",
      "1                5.4084   4.5166   5.5254   6.6659   5.5392   5.2352   7.8130   \n",
      "2               11.7646  15.7579  10.5553  11.7164  11.1714   9.4911  12.5491   \n",
      "3                5.4338   4.7139   5.2118   6.6280   5.4487   5.5255   8.2078   \n",
      "4               11.9036  13.7989  11.7388  11.6912  11.3379   9.1987  15.1357   \n",
      "5                5.4810   4.8976   5.2497   6.6624   5.5921   5.4969   8.1883   \n",
      "6               11.2929  13.7584   9.1466  15.4543  13.0933  12.3646  11.7948   \n",
      "7                8.1014   9.6989   6.3436  10.8129  10.5661   7.6788   7.9391   \n",
      "8               16.1052  18.4056  12.9217  13.6548  24.0948  27.4331  19.9226   \n",
      "9                8.6168   9.5626   7.3891  10.5307  11.9329   7.7953   8.9230   \n",
      "10              16.1038  18.4510  12.9237  13.6858  23.9510  27.5083  19.9355   \n",
      "11               8.6295   9.5959   7.3687  10.5836  11.8951   7.7050   8.9143   \n",
      "12              17.1838  17.4118  14.5725  21.9899  16.4564  16.5139  19.2928   \n",
      "13              13.5982  13.3074  11.8648  15.8414  12.8857  13.5557  15.9971   \n",
      "14              22.8140  24.5202  18.8081  24.1922  25.5158  21.5408  20.1852   \n",
      "15              16.6810  14.1044  14.1341  19.0642  13.5756  21.0341  20.9458   \n",
      "16              22.8186  24.4883  18.8617  23.9648  25.7620  21.5715  20.1280   \n",
      "17              19.5550  20.0743  13.7492  24.0272  17.1300  21.4882  20.0996   \n",
      "18              20.9537  20.9535  21.4743  21.8783  23.8256  22.4091  21.3453   \n",
      "19              17.4903  16.8401  17.5854  18.5029  19.3390  18.4670  18.3508   \n",
      "20              25.1193  28.1063  26.4237  23.5065  21.7771  27.1311  31.5921   \n",
      "21              23.8667  27.3983  26.8236  21.3144  21.0846  26.1433  31.4758   \n",
      "22              24.8331  27.4793  27.0090  21.5464  21.3895  26.5660  31.1365   \n",
      "23              24.7766  27.5407  26.9974  21.3458  21.1357  26.3250  31.4617   \n",
      "24              27.5329  24.6716  28.2918  27.8276  27.0857  24.3790  28.0369   \n",
      "25              23.8712  21.5920  23.7644  25.3051  23.2405  21.7877  23.5777   \n",
      "26              31.6410  28.5178  40.3921  28.9896  30.8735  30.2317  34.4885   \n",
      "27              31.6656  28.6379  39.9663  28.9248  30.8432  30.3934  34.4544   \n",
      "28              31.6409  28.5223  40.4056  28.9862  30.8728  30.2419  34.4910   \n",
      "29              31.6666  28.5939  39.9811  28.9621  30.8600  30.3604  34.4508   \n",
      "30              32.1733  33.4833  31.4847  39.3470  33.9296  33.4922  37.6832   \n",
      "31              28.4237  29.2742  27.6304  31.8280  28.7715  27.7590  31.3330   \n",
      "32              36.7787  35.2172  37.7392  37.1885  37.2329  42.1517  40.1436   \n",
      "33              36.7641  35.2148  37.7634  37.1935  37.1552  42.1112  40.2419   \n",
      "34              36.7772  35.2046  37.7466  37.1705  37.2266  42.1245  40.1198   \n",
      "35              36.7598  35.1640  37.8137  36.9947  37.1064  42.0680  39.9224   \n",
      "36              36.8255  40.2658  38.6375  35.3971  35.5871  36.2119  36.7560   \n",
      "37              33.4744  36.0544  34.3520  34.2692  31.7140  32.1034  34.9014   \n",
      "38              41.2319  37.9423  53.0953  36.6010  44.6575  40.3041  41.4450   \n",
      "39              41.2197  37.9049  52.7816  36.5857  44.5950  40.3400  41.4840   \n",
      "40              41.2222  37.8658  52.7468  36.5892  44.7604  40.2977  41.4568   \n",
      "41              41.1925  37.6738  51.9916  36.4847  45.3160  40.2346  41.5227   \n",
      "42              42.8217  41.0421  42.8148  41.2278  38.9211  46.4598  49.5841   \n",
      "43              39.1507  38.1748  38.3433  38.0858  35.6295  42.1483  42.8090   \n",
      "44              45.3085  44.5019  48.0592  48.4297  45.0597  48.8957  42.5576   \n",
      "45              45.2928  44.5172  47.9029  48.6373  45.0336  48.6752  42.5396   \n",
      "46              45.1389  44.5362  47.7384  48.6874  45.5312  48.2858  42.2310   \n",
      "47              45.0668  44.4746  47.4797  48.5105  45.7581  48.1830  42.3136   \n",
      "48              47.1896  48.3298  46.2606  45.9491  48.2987  47.9214  47.3137   \n",
      "49              43.3193  43.2275  42.8015  43.7083  43.4803  43.9305  43.3441   \n",
      "50              53.9295  54.6444  51.8571  49.3119  56.0809  53.3861  58.2839   \n",
      "51              53.8569  54.3512  51.8486  49.2354  56.1120  53.3770  58.3472   \n",
      "52              53.9046  54.6369  51.8847  49.3142  55.9283  53.3869  58.2342   \n",
      "53              53.7181  54.3721  51.8719  49.0444  55.4161  53.3379  58.0710   \n",
      "\n",
      "          6        7        8        9  \n",
      "0    6.6378  10.3879   9.5944   8.2098  \n",
      "1    4.5757   4.2648   6.3294   4.8290  \n",
      "2   10.4448  13.1749  10.8998  16.9683  \n",
      "3    4.8871   3.9670   6.3999   4.7609  \n",
      "4    9.4022  15.3920  13.0271  12.9398  \n",
      "5    4.8957   4.0179   6.4132   4.7595  \n",
      "6    9.8675   9.9648  11.1463  10.1631  \n",
      "7    6.6366   6.9418   8.2472   8.5982  \n",
      "8   16.4114  16.2622  18.5485  12.1995  \n",
      "9    6.8956   7.6877   8.1271   9.7234  \n",
      "10  16.3833  16.2322  18.4819  12.2022  \n",
      "11   7.0887   7.6888   8.1044   9.6892  \n",
      "12  17.8833  16.7959  16.0282  18.8453  \n",
      "13  13.0811  13.9051  12.9423  14.2603  \n",
      "14  40.5173  22.5952  21.8889  34.5970  \n",
      "15  17.5315  19.1848  14.5868  32.6354  \n",
      "16  40.2060  22.5963  21.9296  34.5048  \n",
      "17  22.1809  22.4939  21.8855  33.6799  \n",
      "18  19.9572  17.9903  22.7495  19.9439  \n",
      "19  15.9029  15.6215  17.7987  18.0180  \n",
      "20  24.6986  21.8712  29.7768  26.4432  \n",
      "21  23.2293  19.3057  33.2092  27.2891  \n",
      "22  23.2695  23.2151  33.7433  27.1200  \n",
      "23  23.2224  23.5311  33.6342  27.3875  \n",
      "24  29.0369  31.5116  27.3450  38.7471  \n",
      "25  24.8055  26.4793  23.0828  31.5581  \n",
      "26  34.8872  35.1768  35.2483  29.8099  \n",
      "27  34.7806  35.2715  35.3608  29.7951  \n",
      "28  34.8866  35.1608  35.2543  29.7996  \n",
      "29  34.8106  35.3126  35.2971  29.8261  \n",
      "30  30.8843  30.1837  31.8495  28.7370  \n",
      "31  28.0110  27.0621  28.6695  26.6408  \n",
      "32  42.2767  35.7544  36.7600  32.5434  \n",
      "33  42.2798  35.7976  36.7149  32.4913  \n",
      "34  42.3219  35.7650  36.7646  32.5463  \n",
      "35  42.5510  35.8565  36.7458  32.5678  \n",
      "36  37.3441  38.9840  39.0899  34.0119  \n",
      "37  33.1515  33.9515  33.8705  32.2647  \n",
      "38  41.9200  58.6020  58.4123  39.8905  \n",
      "39  41.8526  57.3597  57.6692  39.9431  \n",
      "40  41.9160  57.6410  58.1010  39.9689  \n",
      "41  41.8183  54.0544  56.4216  40.3865  \n",
      "42  45.3551  46.0651  40.5241  49.3788  \n",
      "43  41.3196  40.9364  38.0386  41.9715  \n",
      "44  42.6766  49.5480  47.2683  43.4043  \n",
      "45  42.6448  49.3749  47.4701  43.3852  \n",
      "46  42.3355  49.0490  47.6991  42.9897  \n",
      "47  42.1464  48.6514  47.7717  42.8443  \n",
      "48  45.7818  47.8598  46.8519  48.4407  \n",
      "49  42.2400  43.5394  42.9388  44.3609  \n",
      "50  52.2783  58.2617  56.4137  62.4934  \n",
      "51  52.2647  58.1476  56.2102  61.8980  \n",
      "52  52.2368  58.1661  56.3666  62.0855  \n",
      "53  52.2480  57.6784  56.2579  61.1668  \n"
     ]
    }
   ],
   "source": [
    "data = pd.read_csv('../../public_code/datasets/2d-grid/2d-grid-simulation-19-02-12-23-54.csv',\n",
    "                       skiprows=3,\n",
    "                      names=['input_snr (dB)', 'sigma_sq', 'penalty_f', 'gamma_mult', 'rho_mult', 'penalty_param',\n",
    "                              'avg output_snr (dB)'] + [str(e) for e in range(10)])\n",
    "print(data)\n",
    "\n",
    "input_snr = np.array(sorted(list(set(data['input_snr (dB)'].values))))\n",
    "\n",
    "L1 = data.loc[data['penalty_f'] == 'L1', ['gamma_mult', 'rho_mult']].values\n",
    "SCAD = data.loc[data['penalty_f'] == 'SCAD-L1', ['gamma_mult', 'rho_mult']].values\n",
    "MCP = data.loc[data['penalty_f'] == 'MCP-L1', ['gamma_mult', 'rho_mult']].values\n",
    "    \n",
    "L1_vec = L1[range(0,len(L1),2)]\n",
    "SCAD_vec = SCAD[range(0,len(SCAD),2)]\n",
    "MCP_vec = MCP[range(0,len(MCP),2)]\n",
    "    \n",
    "L1 = L1[range(1,len(L1),2)]\n",
    "SCAD = SCAD[range(1,len(SCAD),2)]\n",
    "MCP = MCP[range(1,len(MCP),2)]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "n = 400 m = 760\n",
      "sparsity 24\n"
     ]
    }
   ],
   "source": [
    "n = 20\n",
    "Y_HIGH = 10\n",
    "Y_LOW = -5\n",
    "k = 0\n",
    "Gnx, signal_2d, y_true, xs, ys = create2DSignal(k, n, Y_HIGH=Y_HIGH, Y_LOW=Y_LOW)\n",
    "pspace = (\n",
    "    hp.loguniform('gamma_mult', -4, 4),\n",
    "    hp.loguniform('rho_mult', -3, 3)\n",
    ")\n",
    "max_evals = 30\n",
    "d = 10\n",
    "n = nx.number_of_nodes(Gnx)\n",
    "\n",
    "Dk = penalty_matrix(Gnx, k)\n",
    "print ('n =', n,  'm =', Gnx.number_of_edges())\n",
    "print ('sparsity', np.count_nonzero(Dk*y_true))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "DTD = Dk.T.dot(Dk).toarray()\n",
    "[S, V] = np.linalg.eigh(DTD)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "20\n",
      "sigma squared: 0.11249999999999999\n",
      "input signal MSE: 45.08730299353278\n",
      "input signal SNR: 19.991582560704053\n",
      "L1 vec MSE: 1.1105974454221341\n",
      "L1 vec SNR: 36.076558435946275\n",
      "100%|██████████| 30/30 [00:07<00:00,  1.73it/s, best loss: -38.67422783323206]\n",
      "100%|██████████| 30/30 [00:06<00:00,  2.03it/s, best loss: -38.66861109609185] \n",
      "sigma squared: 0.11249999999999999\n",
      "input signal MSE: 43.8320595837574\n",
      "input signal SNR: 20.114206359836725\n",
      "L1 vec MSE: 0.8679640278729892\n",
      "L1 vec SNR: 37.14710787239393\n",
      "100%|██████████| 30/30 [00:08<00:00,  1.07it/s, best loss: -40.210219737928654]\n",
      "100%|██████████| 30/30 [00:04<00:00,  6.81it/s, best loss: -39.96758841863024] \n",
      "sigma squared: 0.11249999999999999\n",
      "input signal MSE: 45.267429419618416\n",
      "input signal SNR: 19.97426680657877\n",
      "L1 vec MSE: 0.8210842087069224\n",
      "L1 vec SNR: 37.38824814023335\n",
      "100%|██████████| 30/30 [00:07<00:00,  3.80it/s, best loss: -41.54619294417688]\n",
      "100%|██████████| 30/30 [00:04<00:00,  5.17it/s, best loss: -41.36650236536152] \n",
      "sigma squared: 0.11249999999999999\n",
      "input signal MSE: 45.57223351455486\n",
      "input signal SNR: 19.945121997139413\n",
      "L1 vec MSE: 0.8631131637891141\n",
      "L1 vec SNR: 37.17144773463146\n",
      "100%|██████████| 30/30 [00:06<00:00,  1.83it/s, best loss: -40.457430846848624]\n",
      "100%|██████████| 30/30 [00:05<00:00,  5.11it/s, best loss: -40.43376226592328] \n",
      "sigma squared: 0.11249999999999999\n",
      "input signal MSE: 43.326095189321144\n",
      "input signal SNR: 20.164629642320374\n",
      "L1 vec MSE: 1.0608325118344255\n",
      "L1 vec SNR: 36.27565692482822\n",
      "100%|██████████| 30/30 [00:04<00:00,  8.14it/s, best loss: -42.13748671135502]\n",
      "100%|██████████| 30/30 [00:06<00:00,  7.71it/s, best loss: -41.428605213590785]\n",
      "sigma squared: 0.11249999999999999\n",
      "input signal MSE: 46.15068646447885\n",
      "input signal SNR: 19.89034348488864\n",
      "L1 vec MSE: 1.0138009054782748\n",
      "L1 vec SNR: 36.472598389968155\n",
      "100%|██████████| 30/30 [00:03<00:00,  8.67it/s, best loss: -42.196957318024566]\n",
      "100%|██████████| 30/30 [00:04<00:00,  4.39it/s, best loss: -42.19694291750616]\n",
      "sigma squared: 0.11249999999999999\n",
      "input signal MSE: 44.11522050565744\n",
      "input signal SNR: 20.086240594282586\n",
      "L1 vec MSE: 1.1034357914019342\n",
      "L1 vec SNR: 36.104654470081975\n",
      "100%|██████████| 30/30 [00:05<00:00,  5.28it/s, best loss: -41.4786019496933] \n",
      "100%|██████████| 30/30 [00:08<00:00,  3.27it/s, best loss: -41.475732218545076]\n",
      "sigma squared: 0.11249999999999999\n",
      "input signal MSE: 45.22089951232674\n",
      "input signal SNR: 19.978733169036733\n",
      "L1 vec MSE: 0.8648903654304767\n",
      "L1 vec SNR: 37.16251454525821\n",
      "100%|██████████| 30/30 [00:05<00:00,  5.39it/s, best loss: -40.34114351538311]\n",
      "100%|██████████| 30/30 [00:05<00:00,  3.73it/s, best loss: -40.2331936063429] \n",
      "sigma squared: 0.11249999999999999\n",
      "input signal MSE: 45.43268417095266\n",
      "input signal SNR: 19.95844118008548\n",
      "L1 vec MSE: 1.100307435256242\n",
      "L1 vec SNR: 36.11698466089266\n",
      "100%|██████████| 30/30 [00:05<00:00,  4.21it/s, best loss: -40.529626443598005]\n",
      "100%|██████████| 30/30 [00:03<00:00,  7.21it/s, best loss: -40.54279900095245] \n",
      "sigma squared: 0.11249999999999999\n",
      "input signal MSE: 45.74824577752749\n",
      "input signal SNR: 19.928380681323578\n",
      "L1 vec MSE: 1.1737698752184644\n",
      "L1 vec SNR: 35.83629554614381\n",
      "100%|██████████| 30/30 [00:04<00:00,  2.92it/s, best loss: -40.06107847648291]\n",
      "100%|██████████| 30/30 [00:03<00:00,  8.93it/s, best loss: -40.003445383145966]\n"
     ]
    }
   ],
   "source": [
    "num_trial = 10\n",
    "idx = 6\n",
    "\n",
    "INPUT_SNR = input_snr[idx]\n",
    "print(INPUT_SNR)\n",
    "\n",
    "Ys = []\n",
    "scad_opts = []\n",
    "mcp_opts = []\n",
    "\n",
    "y_true_norm_sq = np.linalg.norm(y_true, 2) ** 2\n",
    "sigma_sq = y_true_norm_sq/(10**(INPUT_SNR/10.0))/n\n",
    "\n",
    "for trial in range(num_trial):\n",
    "    Y = np.tile(y_true.reshape(-1, 1), (1, d)) + np.random.normal(scale=np.sqrt(sigma_sq), size=(n, d))\n",
    "    Ys.append(Y)\n",
    "    ses = np.linalg.norm(Y - np.tile(y_true.reshape(-1, 1), (1, d)), 2, axis=0) ** 2\n",
    "\n",
    "    mse = np.mean(ses)\n",
    "    mean_snr = 10*np.log10(y_true_norm_sq/mse)\n",
    "    \n",
    "    print ('sigma squared:', sigma_sq)\n",
    "    print ('input signal MSE:', mse)\n",
    "    print ('input signal SNR:', mean_snr)\n",
    "\n",
    "# opt_param, B_hat, snr_avg, penalty_param, output_snr = autotune_denoising(Gnx, Y, y_true.reshape(-1,1), 0,\n",
    "#                                                                                           'L1', sigma_sq,\n",
    "#                                                                                           max_evals, pspace,\n",
    "#                                                                                           B_init=None, vec=True)\n",
    "\n",
    "    res = admm_denoising_snr(L1_vec[idx], Y, sigma_sq, y_true, Dk, 'L1', eig=(S, V), B_init=None, vec=True)\n",
    "\n",
    "    B_hat = res['B']\n",
    "    l1_snr = -res['loss']\n",
    "    mse = res['mse']\n",
    "    \n",
    "    print ('L1 vec MSE:', mse)\n",
    "    print ('L1 vec SNR:', l1_snr)\n",
    "\n",
    "    l1_B = B_hat.copy()\n",
    "\n",
    "    opt_param, B_hat, snr_avg, penalty_param, output_snr = autotune_denoising(Y, y_true.reshape(-1,1), Dk,\n",
    "                                                                                          'SCAD', sigma_sq,\n",
    "                                                                                          max_evals, pspace,eig=(S, V),\n",
    "                                                                                          B_init=l1_B, vec=True)\n",
    "    scad_opts.append(opt_param.copy())\n",
    "     \n",
    "    opt_param, B_hat, snr_avg, penalty_param, output_snr = autotune_denoising(Y, y_true.reshape(-1,1), Dk,\n",
    "                                                                                          'MCP', sigma_sq,\n",
    "                                                                                          max_evals, pspace,eig=(S, V),\n",
    "                                                                                          B_init=l1_B, vec=True)\n",
    "    mcp_opts.append(opt_param.copy())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[{'gamma_mult': 51.96405443765451, 'rho_mult': 15.336395908647924},\n",
       " {'gamma_mult': 54.57974044599013, 'rho_mult': 7.44397108703288},\n",
       " {'gamma_mult': 52.26544051099, 'rho_mult': 4.116563823124163},\n",
       " {'gamma_mult': 28.79629752761651, 'rho_mult': 2.9485304776220573},\n",
       " {'gamma_mult': 38.51848553970784, 'rho_mult': 2.0855839079096885},\n",
       " {'gamma_mult': 47.38145912222227, 'rho_mult': 1.6037157011011094},\n",
       " {'gamma_mult': 11.30821498549191, 'rho_mult': 9.545474045827463},\n",
       " {'gamma_mult': 53.33023916010095, 'rho_mult': 6.545485387620458},\n",
       " {'gamma_mult': 51.20428054765246, 'rho_mult': 1.4418317593845495},\n",
       " {'gamma_mult': 39.598618257085896, 'rho_mult': 4.928404028071851}]"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "[{'gamma_mult': 53.72146631645505, 'rho_mult': 3.9249656869531298},\n",
       " {'gamma_mult': 22.83418675024217, 'rho_mult': 1.7393301782014685},\n",
       " {'gamma_mult': 52.87736306025207, 'rho_mult': 1.3208463855453496},\n",
       " {'gamma_mult': 23.036806304373204, 'rho_mult': 2.5575559027538395},\n",
       " {'gamma_mult': 54.33436802851751, 'rho_mult': 0.4307624510344473},\n",
       " {'gamma_mult': 23.41707953931136, 'rho_mult': 3.2767788397804276},\n",
       " {'gamma_mult': 25.068518797650302, 'rho_mult': 3.697252727689832},\n",
       " {'gamma_mult': 50.83713503349903, 'rho_mult': 2.6930947562576875},\n",
       " {'gamma_mult': 52.82053988049464, 'rho_mult': 1.845052626614221},\n",
       " {'gamma_mult': 39.9743018632703, 'rho_mult': 1.2133909448276117}]"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "sigma squared: 0.11249999999999999\n",
      "input signal MSE: 45.74824577752749\n",
      "input signal SNR: 19.928380681323578\n",
      "L1 0.7852275 11.470289228999999 0\n",
      "L1\n",
      "output signal MSE: 1.1105974454221341\n",
      "output signal SNR: 36.076558435946275\n",
      "iteration: (102,)\n",
      "SCAD 5.845956124236132 89.65589758587029 3.7\n",
      "SCAD\n",
      "output signal MSE: 0.6106463006555811\n",
      "output signal SNR: 38.67422783323206\n",
      "iteration: (108,)\n",
      "MCP 6.043664960601193 23.721177593800622 1.4\n",
      "MCP\n",
      "output signal MSE: 0.6114365613575641\n",
      "output signal SNR: 38.66861109609185\n",
      "iteration: (32,)\n",
      "L1 0.7852275 11.470289228999999 0\n",
      "L1\n",
      "output signal MSE: 0.8679640278729892\n",
      "output signal SNR: 37.14710787239393\n",
      "iteration: (114,)\n",
      "SCAD 6.140220800173889 45.70762610449233 3.7\n",
      "SCAD\n",
      "output signal MSE: 0.4287365806794883\n",
      "output signal SNR: 40.210219737928654\n",
      "iteration: (61,)\n",
      "MCP 2.568846009402244 4.468071387305736 1.4\n",
      "MCP\n",
      "output signal MSE: 0.45337093214349417\n",
      "output signal SNR: 39.96758841863024\n",
      "iteration: (17,)\n",
      "L1 0.7852275 11.470289228999999 0\n",
      "L1\n",
      "output signal MSE: 0.8210842087069224\n",
      "output signal SNR: 37.38824814023335\n",
      "iteration: (104,)\n",
      "SCAD 5.879862057486374 24.204827430808812 3.7\n",
      "SCAD\n",
      "output signal MSE: 0.3152050881145834\n",
      "output signal SNR: 41.546192944176866\n",
      "iteration: (34,)\n",
      "MCP 5.948703344278358 7.857323310971602 1.4\n",
      "MCP\n",
      "output signal MSE: 0.3285203502978173\n",
      "output signal SNR: 41.36650236536152\n",
      "iteration: (20,)\n",
      "L1 0.7852275 11.470289228999999 0\n",
      "L1\n",
      "output signal MSE: 0.8631131637891141\n",
      "output signal SNR: 37.17144773463146\n",
      "iteration: (110,)\n",
      "SCAD 3.239583471856857 9.55201060157062 3.7\n",
      "SCAD\n",
      "output signal MSE: 0.40501343437854703\n",
      "output signal SNR: 40.45743084684862\n",
      "iteration: (19,)\n",
      "MCP 2.591640709241985 6.628265993738986 1.4\n",
      "MCP\n",
      "output signal MSE: 0.4072267395788668\n",
      "output signal SNR: 40.43376226592328\n",
      "iteration: (18,)\n",
      "L1 0.7852275 11.470289228999999 0\n",
      "L1\n",
      "output signal MSE: 1.0608325118344255\n",
      "output signal SNR: 36.27565692482822\n",
      "iteration: (98,)\n",
      "SCAD 4.333329623217131 9.037522529850003 3.7\n",
      "SCAD\n",
      "output signal MSE: 0.2750830573960819\n",
      "output signal SNR: 42.13748671135501\n",
      "iteration: (20,)\n",
      "MCP 6.112616403208219 2.63308562407934 1.4\n",
      "MCP\n",
      "output signal MSE: 0.32385603347130126\n",
      "output signal SNR: 41.428605213590785\n",
      "iteration: (13,)\n",
      "L1 0.7852275 11.470289228999999 0\n",
      "L1\n",
      "output signal MSE: 1.0138009054782748\n",
      "output signal SNR: 36.472598389968155\n",
      "iteration: (99,)\n",
      "SCAD 5.330414151250005 8.548468867731176 3.7\n",
      "SCAD\n",
      "output signal MSE: 0.2713418501827192\n",
      "output signal SNR: 42.196957318024566\n",
      "iteration: (19,)\n",
      "MCP 2.6344214481725277 8.632416456435449 1.4\n",
      "MCP\n",
      "output signal MSE: 0.2713427499108873\n",
      "output signal SNR: 42.19694291750616\n",
      "iteration: (19,)\n",
      "L1 0.7852275 11.470289228999999 0\n",
      "L1\n",
      "output signal MSE: 1.1034357914019342\n",
      "output signal SNR: 36.104654470081975\n",
      "iteration: (102,)\n",
      "SCAD 1.2721741858678397 12.143505672973147 3.7\n",
      "SCAD\n",
      "output signal MSE: 0.320149124716213\n",
      "output signal SNR: 41.4786019496933\n",
      "iteration: (20,)\n",
      "MCP 2.8202083647356586 10.427023069172595 1.4\n",
      "MCP\n",
      "output signal MSE: 0.3203607427688457\n",
      "output signal SNR: 41.47573221854507\n",
      "iteration: (19,)\n",
      "L1 0.7852275 11.470289228999999 0\n",
      "L1\n",
      "output signal MSE: 0.8648903654304767\n",
      "output signal SNR: 37.16251454525821\n",
      "iteration: (105,)\n",
      "SCAD 5.999651905511357 39.27063387833382 3.7\n",
      "SCAD\n",
      "output signal MSE: 0.4160046281063773\n",
      "output signal SNR: 40.3411435153831\n",
      "iteration: (53,)\n",
      "MCP 5.719177691268641 15.402287450461523 1.4\n",
      "MCP\n",
      "output signal MSE: 0.4264745829175112\n",
      "output signal SNR: 40.2331936063429\n",
      "iteration: (27,)\n",
      "L1 0.7852275 11.470289228999999 0\n",
      "L1\n",
      "output signal MSE: 1.100307435256242\n",
      "output signal SNR: 36.11698466089266\n",
      "iteration: (92,)\n",
      "SCAD 5.760481561610901 8.305645264879704 3.7\n",
      "SCAD\n",
      "output signal MSE: 0.3983362856628539\n",
      "output signal SNR: 40.529626443598005\n",
      "iteration: (20,)\n",
      "MCP 5.942310736555647 10.963876032639883 1.4\n",
      "MCP\n",
      "output signal MSE: 0.3971299249230774\n",
      "output signal SNR: 40.54279900095245\n",
      "iteration: (25,)\n",
      "L1 0.7852275 11.470289228999999 0\n",
      "L1\n",
      "output signal MSE: 1.1737698752184644\n",
      "output signal SNR: 35.83629554614381\n",
      "iteration: (107,)\n",
      "SCAD 4.454844553922163 21.955273843983935 3.7\n",
      "SCAD\n",
      "output signal MSE: 0.4437155676802826\n",
      "output signal SNR: 40.06107847648291\n",
      "iteration: (33,)\n",
      "MCP 4.497108959617909 5.456751289503492 1.4\n",
      "MCP\n",
      "output signal MSE: 0.4496431436167172\n",
      "output signal SNR: 40.003445383145966\n",
      "iteration: (18,)\n"
     ]
    }
   ],
   "source": [
    "scad_opts\n",
    "mcp_opts\n",
    "\n",
    "idx = 6 # Input SNR = 20dB\n",
    "\n",
    "l1_times = []\n",
    "scad_times = []\n",
    "mcp_times = []\n",
    "l1_snrs = []\n",
    "scad_snrs = []\n",
    "mcp_snrs = []\n",
    "\n",
    "ses = np.linalg.norm(Y - np.tile(y_true.reshape(-1, 1), (1, d)), 2, axis=0) ** 2\n",
    "\n",
    "mse = np.mean(ses)\n",
    "mean_snr = 10*np.log10(y_true_norm_sq/mse)\n",
    "    \n",
    "print ('sigma squared:', sigma_sq)\n",
    "print ('input signal MSE:', mse)\n",
    "print ('input signal SNR:', mean_snr)\n",
    "\n",
    "for trial in range(num_trial):\n",
    "    scad_opt = scad_opts[trial]\n",
    "    mcp_opt = mcp_opts[trial]\n",
    "    Y = Ys[trial]\n",
    "    PARAM_GRID = [\n",
    "            {'penalty_f': [\"L1\"], 'gamma_ratio': [L1_vec[idx,0]], 'rho_ratio': [L1_vec[idx,1]], 'penalty_param': [0]},\n",
    "            {'penalty_f': [\"SCAD\"], 'gamma_ratio': [scad_opt['gamma_mult']], 'rho_ratio': [scad_opt['rho_mult']], 'penalty_param': [3.7]},\n",
    "            {'penalty_f': [\"MCP\"], 'gamma_ratio': [mcp_opt['gamma_mult']], 'rho_ratio': [mcp_opt['rho_mult']], 'penalty_param': [1.4]}\n",
    "        ]\n",
    "    for params in ParameterGrid(PARAM_GRID):\n",
    "        penalty_f = params['penalty_f']\n",
    "        gamma = params['gamma_ratio']*sigma_sq\n",
    "        penalty_param = params['penalty_param']\n",
    "\n",
    "        if penalty_f == 'L1':\n",
    "            l1_beta = None\n",
    "            rho = params['rho_ratio']*gamma\n",
    "        else:\n",
    "            rho = max(params['rho_ratio']*gamma, 1/penalty_param)\n",
    "\n",
    "        print (penalty_f, gamma, rho, penalty_param)\n",
    "        \n",
    "        t = time()\n",
    "        invX = V.dot(np.diag(1/(1+rho*S))).dot(V.T)\n",
    "        B, obj, err_path = admm(Y, gamma, rho, Dk, penalty_f, penalty_param,max_iter=500, tol_abs=10 ** (-3), \n",
    "                               tol_rel=10 ** (-2), B_init=l1_beta, invX=invX)\n",
    "        t1 = time()-t\n",
    "        \n",
    "        ses = np.linalg.norm(B - np.tile(y_true.reshape(-1, 1), (1, d)), 2, axis=0) ** 2\n",
    "        mse = np.mean(ses)\n",
    "        snr = 10*np.log10(y_true_norm_sq/mse)\n",
    "\n",
    "        if penalty_f == 'L1':\n",
    "            l1_times.append(t1)\n",
    "            l1_beta = B\n",
    "            l1_snrs.append(snr)\n",
    "        elif penalty_f == 'SCAD':\n",
    "            scad_times.append(t1)\n",
    "            scad_snrs.append(snr)\n",
    "        else:\n",
    "            mcp_times.append(t1)\n",
    "            mcp_snrs.append(snr)            \n",
    "\n",
    "        print (penalty_f)\n",
    "        print ('output signal MSE:', mse)\n",
    "        print ('output signal SNR:', snr)\n",
    "        print ('iteration:', np.shape(err_path[3]))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[0.3338279724121094,\n",
       " 0.37367796897888184,\n",
       " 0.2660641670227051,\n",
       " 0.4643402099609375,\n",
       " 0.2235116958618164,\n",
       " 0.29561710357666016,\n",
       " 0.22907185554504395,\n",
       " 0.5302920341491699,\n",
       " 0.3063070774078369,\n",
       " 0.25487208366394043]"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "[0.32982897758483887,\n",
       " 0.2277512550354004,\n",
       " 0.13531112670898438,\n",
       " 0.08294796943664551,\n",
       " 0.14155197143554688,\n",
       " 0.09302902221679688,\n",
       " 0.1179499626159668,\n",
       " 0.1589641571044922,\n",
       " 0.1616981029510498,\n",
       " 0.2210710048675537]"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "[0.12143301963806152,\n",
       " 0.0743408203125,\n",
       " 0.07945013046264648,\n",
       " 0.07322001457214355,\n",
       " 0.19911503791809082,\n",
       " 0.07359623908996582,\n",
       " 0.08317899703979492,\n",
       " 0.10697579383850098,\n",
       " 0.09918713569641113,\n",
       " 0.10935091972351074]"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "[36.076558435946275,\n",
       " 37.14710787239393,\n",
       " 37.38824814023335,\n",
       " 37.17144773463146,\n",
       " 36.27565692482822,\n",
       " 36.472598389968155,\n",
       " 36.104654470081975,\n",
       " 37.16251454525821,\n",
       " 36.11698466089266,\n",
       " 35.83629554614381]"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "[38.67422783323206,\n",
       " 40.210219737928654,\n",
       " 41.546192944176866,\n",
       " 40.45743084684862,\n",
       " 42.13748671135501,\n",
       " 42.196957318024566,\n",
       " 41.4786019496933,\n",
       " 40.3411435153831,\n",
       " 40.529626443598005,\n",
       " 40.06107847648291]"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "[38.66861109609185,\n",
       " 39.96758841863024,\n",
       " 41.36650236536152,\n",
       " 40.43376226592328,\n",
       " 41.428605213590785,\n",
       " 42.19694291750616,\n",
       " 41.47573221854507,\n",
       " 40.2331936063429,\n",
       " 40.54279900095245,\n",
       " 40.003445383145966]"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "l1_times\n",
    "scad_times\n",
    "mcp_times\n",
    "l1_snrs\n",
    "scad_snrs\n",
    "mcp_snrs\n",
    "np.savez('2d-grid-runtime'+'.npz', l1_times=l1_times, scad_times=scad_times, mcp_times=mcp_times, \n",
    "         l1_snrs=l1_snrs, scad_snrs=scad_snrs, mcp_snrs=mcp_snrs, scad_opts=scad_opts, mcp_opts=mcp_opts, \n",
    "        Ys=Ys, input_snr=INPUT_SNR)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.collections.PathCollection at 0x11ed201d0>"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.collections.PathCollection at 0x11ed20748>"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.collections.PathCollection at 0x11ed205c0>"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "Text(0.5, 0, 'Gain in SNR (dB)')"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "Text(0, 0.5, 'Runtime (sec)')"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.legend.Legend at 0x11ed20f28>"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAVAAAAFGCAYAAAA8Z8dfAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4xLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvDW2N/gAAIABJREFUeJzt3X1wHPd5H/DvcwBJARSkAynQYmWn1EFySllIQwBShw5rKdFRneFInakC6OR0Ms4UIaB2EttxE0LOVAyn8lQGnRlXo7gyILEZ0tO6IGBnFMlKbYKJJcqELQJHMZiKjkmcXsyhZNIgjwZJEG/39I/bPewd9l6wuNvdu/t+ZjDA7S7ufgB5D35v+zyiqiAiopULeN0AIqJyxQBKROQQAygRkUMMoEREDjGAEhE5xABKROQQAygRkUMMoEREDjGAEhE5VOt1A1bjtttu0y1btnjdDCKqMOPj479U1aZ815V1AN2yZQvGxsa8bgYRVRgReb+Q6ziEJyJyiAGUiMghBlAiIocYQImIHCrrRaRs5ufnce7cOdy4ccPrppBP3XTTTfj4xz+ONWvWeN0UKmOuB1AR6QAQB9CqqvszzrUCGAcQMw6NqGrPSl/j3LlzaGhowJYtWyAiq24zVRZVxdTUFM6dO4c777zT6+ZQGXM1gBoBEqo6IiIhEWlV1ajlkg2qKpZr405e58aNGwyelJWIYOPGjbh48aLXTaEy5/YcaARLQTEGIGw9qaojloftqhqDQwyelAv/f1AxuB1AgwAuWR5vtLtIRMIADmc51y0iYyIyxh4EVbp8NctY08xbfl2F36mqtsN3VR1Q1XZVbW9qynunFVHZmj15ALNvPZc1SKoqZt96DrMnD7jcMjK5HUDjADYYXwcBTGW5rtWd5lSvaDSKnTt3orm5Gc3Nzejt7QUA7N+/P3WssbERjY2Nqcf79yfX/KzHzI94vPDp6sbGxqzXj4yMpF6nmqkqdG4a86eHbIOoGTznTw9B56bZE/WI26vwgwDaja9DAEYAQESCZo9TREIut6lsDAwMoLu7e9XPE41G0dnZiaGhIbS2Jv9WDQ8PAwD27NmDPXv2AEAqqPb19S17jvHxcQSDwZyvMzIygnA4nPMaq87OTsRisRV9T6USEay7/wsAgPnTQwCAdfd/ASKSFjzXbO1MHSf3udoDNVfcjTnOuGUF/mjGpY4Xj1ZremEBL50/j97JSbx0/jymFxa8asoy/f39RXme3bt3o6+vLxU8AaCjo6Moz20VjUbzX2QxNDSEnp4V71qrWGYQXbO1M60nyuDpH67PgRpzmCOqOmA51mb5OuZk72cxvBmP447RUXzx7Fns//nP8cWzZ3HH6CjeXMHwNJfOzk6MjCxtNGhubgawNGzeuXNnamg7MDCQNrzu6elBNBpFW1tbqmc4MDCAtrY2tLW1pZ43Ho+nruns7FzWhng8jmg0WpKAScWXGUSvHtzB4OkjFXknkhPTCwvYNTGB6cXF1LFriQQAYNfEBM5v346ba1f364pEIhgaGkI4HEY0GkVrayui0ShOnDiByclJxGKxVLDs7+/H5OQkACAWiyEUCmFkZATj4+MAkr27oaGh1OO2tjYcPXo0dS4SidgOvS9duoRQaPWzJA899FDq63A4bPtaVBxmEDWH8gAYPH2CAdQweOECElkm4hOqGLx4EV2bN6/qNTo6OlK9x8HBQUQiEQwODiIWi2Hnzp1LbTHOmewC3uDgYNpwNxKJ4PDhw3j88ccRDAZT85iZNmzYgFhsaYZkeHgYzz77LGKxGC5fvlzwz3L06NFlc6CxWCw1lwoAR44cSTvf0dFRlOBdbcxhu9XsW88xiPoAA6jhzMxMqseZ6VoigbPXrxfldcxe58jICPr6+nDixAl8+ctfThtSm0E2l8xV7KmpqVRA27BhQ+q4uTADJOcYQ6EQWltbUws8HR0d6OjoQGNj46p/tlAotCxwZwvkVBi7OU/zMcCeqNf8ug/UdXfX1WF9wP7XsT4QwF319UV5nUgkgv7+/lSQMx+bzOH34OBg2jEACAaDiMfjiMfj6OzsTH1fPB7H8PCw7eq1OcwfHx9P9f5efPFF9PT0pPVEyX+yLRjZLSyRN9gDNUQ2bcKXjDnHTAERRIq0ab+joyO1hQhI9kg7OzvR1pZcR+vp6UF3dzd6enpSi0w9PT1obW1FOBxGW1sbwuEw+vv7EY1G0dzcjGAwiL6+PoRCoYL2Y7a2tqZWvM0g+vjjjxfl5yuU+fMCyZ7rkSNH0NnZiWg0ikuXLiEajS6bAqgmuVbbc21xIndJOf/1am9vV7uaSKdPn8bWrVtX/HxvxuPYNTGBhCquJRJYHwggIILXWlqwI8+eR1pupftA3eb0/4lbZk8egM5NZw2OZpCVtQ1Yt63LgxZWLhEZV9X2fNexB2qxIxjE+e3bMXjxIs5ev4676usRaWpa9ep7tfJz8CwH67Z1QVWz9izNnih7nt5hZMhwc23tqlfbqThyBY9CzleCfD9fpf/8fsdFJPKlxPWL0GsXcibS0GsXkLjOjFzkHQZQ8h1VBRIJ6I3LtkHUDJ5643LyujKex6fyxiE8+Y6IAOs3AUAySALA+k2pRBpm8JSbGiHGcSIvMICSL2ULogye5CcMoORbmUHUDKQMnuQXnAMlXxMRiBFEU8cYPMknGEDJ18w5z7RjOVbnidzEAGqotuJdXpb0GB4eTn1fW1vbssTLZt7T5uZm3NV8J/Z/7WuQmxoR2PjrkJsa8fAj/w69X/rjZf8mjY2NaGtrw86dO9HW1sbkzFR6qlq2H21tbWrnnXfesT2ezY3oSzrz469rIpGwPZ9IJHTmx1/XG9GXVvS8xdbf31+U5xkfH9dQKKTj4+OpY0NDQ8uu27Nnj+7Zs2fZ8WAwqJcvX877OkeOHFl2bHJyMu37JycndXJyclnbzp49q4vTH+nCxdN64tgPUv82Z8+e1a/+1/+irb9xjy5Of5T2b5bZrv7+fm1tbc3avpX+P6HqAWBMC4hBVd8D1TIq3lUJJT3Mmkdm6r1QKJSWI3T37t345je/iTs/dnNqtb3tt8KpOc+BgQF0/N4foK2tDeM/eTPncL67uxsbNmxIqwJAVExVH0BzpQezBs9ilFBgSQ+gvb0dIyMj6O3ttU2nZ04tIBCwXW2PxWJobm5G52d/Hy/+r5eT1+X4N7FmviIqtqoPoIB7xbvMkh4AbEt69Pf3o7e3F9FoNFXSY3JyMlXiIxQKYXx8HH19fWklPcbHx9Hb25sKvtFoFBs3brQNHMUs6WEG70ISQJuCwSDGx8dTgd76R8MsXQIAgfom2+Bpnt+5cyeOvv4jBOpzpxkMhULMe0olw32ghswci2aexWIW72JJj47UkL2/vz/1B6O3tzf1B8LatszfeX9/P0ZGRlK/KzNvqHUqIpM16BIVGwOohRvFu1jSI11PT0/aVIO1bZmi0WiqiB6QDP5mIM7myJEjaX+MysmVG1fw6QOfxvGu47j1plu9bg7Z4BDewhy2WxW7ZEK1l/QYGRnBwECqojX6+vrS2m3XtuHhYdueZDgcxuHDh7O+1sDAAGKxWNmWcH71Z6/inV++g++d+Z7XTaEsGEANmXOeN3/uzZLUneno6MDAwEBq+G0t6dHW1oaxsTG0tramSno0NzenFojMkh69vb0Ih8OpfZwPPfRQqqRHIawlPczXcKukR3t7O8bHx1N7TC9dupRWEtlsW2dnZ6ptsVgM/f39adMcQPIPyoYNG9JW+x966KHUPlDzD0e5OnjqYPLz2wc9bgllw5IeyL7aXoqFpGrCkh4r893T38UP3/th6vHA+ABmF2exrmYdutu6U8cf3PIgHtv6mActrB4s6VGgXEGSxbtWx8/B04/mF+fxwtgLWEgspB2fXZzF8289DwCoDdRixyd2eNE8slH1AVREIGsbsvYwrUFU1jYweFLJRO6NoOVjLXj024/iw+kPMbMwkzpXV1uHzQ2b8cpnX8E9Tfd42EqyqvoACrB4F/nHPU33YLx7HLftvy3t+NziHKLdUa7G+wwXkQws3kV+cez9Y6hfU4/aQC1qpAa1gVrUr6nHsQ+Oed00ysAASuQzh04dwtW5q9h2+zYc7zqObbdvw9W5qzh06pDXTaMMDKBEPnPm0hnsfWAvRrtGcf8d92O0axR7H9iLM1NnvG4aZeAcKJHPvP3k22mPawI12PfgPux7cJ83DaKs2AMlInKIAZSIyCEG0CrlZUmP5uZm21ylvb29y3Y7WMt7ZLbBWr5jJSn1iIqmkLT1fv0oVkmPTPGZuN7zV/dofCa+qucptkoo6aGqGgqFbEtttLa2ajAYXNbOzJIfdm3I1tZcWNKDsgFLejjn1yw4lVDSwxQOh9Nyh8ZiMbS3p996vHv37lSeUFO23J+RSCTn6xGVAgOojVJlwWFJjyWZKfv6+/uXtTcajRZ8P/2zzz67LFsTUalxGxOWZ8F54/03AACvv/86Pv93n08dX20WHLOkRzgcti3pEYvFUsHSLOkBLGVVHxkZSaVns5b0AIC2tjYcPXo0dS4SiaSliTMVs6SHKRwO275WLpnZ580E06ZCMslb2xCJRPImbyYqNgZQuJcFhyU9OtJ+FnMY39rauuxnLKSWkV0bypnmyMdQyHlyHwMo3M2Cw5IeS3p6elLVOa1/DEy5yntUmtmTB6Bz01mT1qiRdlHWNmDdti4PWkh2OAdqMLPgzC3OpR03s+AUK4VYtZf0sDJ7mUeOHLFte7byHpVGVaFz01mrH5jBc/70EHRuuqglZmh1XO+BikgHgDiAVlXdb3O+FUAIAFTV1XeLmQVnZmEmNVyqq63DsQ+O4ZFPPlKU1+jo6EirVW4t6QEke2Xd3d2pchvmsdbW1lRJj3A4jP7+fkSjUTQ3NyMYDKZKehSyH9Na0sMMTm6V9MgUDoezBnJreQ/z57LrqZa7XIm7rcGTVRF8qJC9TsX6ANAKoMP4uhvJIJp5zZDxeY/deetHsfeBdgx2qOwTvW/gPv3JuZ/ofQP3qewT7Tzc6ej5ql22faB+4bd9oIlEQmd+/HX91V9/Wmd+/HXbx+QOFLgP1O0eaASAubIQAxAGkNq8Z/ROTwCA2vROS83MgvP0Z55GTaAGo12jeOaNZ/DyT192uykVoRrmLospsydq9kbZ8/QvtwNoEMAly+ONGefvA1LD+LBdEBWRbiR7r/i1X/u1ojaOWXDIa2YQNYMnwDpcfubHRaQpVY0CqR5pGlUdUNV2VW1vampyv3VEJaTGnKdVMctqU3G5HUDjAMw9NkEAUxnnp5Ac2pvX3udSu4g8ZwZPc8Ho5s+9iTVbO7OuzpP33B7CDwIwb3gOARgBABEJqmocwDAAs9cZhDEf6oRy0zHl4LdglBk8zWE7y2r7m6s9UMvQPAwgbj4GcNQ4HwMQN4buG9XhNqaamhrMz88Xo8lUoebn51Fb64/7SLIFT2BpTpQ9UX9y/X+Qqg7YHGuzOe94D2gwGMQvfvEL3HHHHQgE/DjNS15KJBL4xS9+gVtv9UeJYBGBrG3Iutpu7YnK2gb2QH1EyvmvWXt7u46NjS07nkgkcO7cOVy7ds2DVlE5WL9+PT7+8Y/76g9svmknTku5R0TGVbU933X+GMMUWSAQKPoWJ6JSyxccGTz9xz9/fomIygwDKBGRQwygREQOMYASETnEAEpE5BADKBGRQwygREQOMYASETnEAEpE5BADKBGRQwygREQOMYASETnEAEpE5BADKFGZyZeCspxTVJYbBlCiMjJ78kDOrPRmdvvZkwdcbll1YgAlKhOqCp2bzlraw1oaROem2RN1QUUmVCaqRLmKzOWqq0SlwwBKVEayBVEGT28wgBKVmcwgagZSBk/3cQ6UqAyJCG78xh/gvvMncSWxAIA1473AAEpUhlQVf/N//xj/ND+D789cBgDWjPcAAyhRmTEXjL71s1cAAIdvuRNrtnZmXZ2n0uEcKFGZ+O7p7+KH7/0Qix+9jcXLZ3B87hoA4PX3X8ee234di1qDxR/14bc/OoXH/+3/5HDeBQygRGVibmEOL5z4BhY0kXZ8dnEWz7/1PACgVgLYfu44Zt96jnOiLlhRABWR3wGwE0ArgA3G4UsAYgDGARxW1V8VtYVka3phAYMXLuDMzAzurqtDZNMmNNTy72Ele6LlCfyLy2fw2I+fw0fz1zGzMJM6V1dbh80Nm/G3T/wtmt89AlnbwODpAilkvkRE/hDAU0gGzREkA+YUgDiAZgBBAO0AtgEYBrBHVd8vUZtT2tvbdWxsrNQv4ztvxuPYNTGBhCquJRJYHwggIILXWlqwIxj0unlUYpdnLqPpa01Y1MXUsRqpwdSeKdx6061QVQbPVRKRcVVtz3dd3i6LiBwGcCuATlU9mefaWwFEAIyISLeq/kOhDabCTC8sYNfEBKYXl9481xLJId2uiQmc374dN7MnWtHe/OBN1K+px8zCTCpY1tXW4dgHx/DIJx9h8HRRzlV4EdkNYFBV/02+4AkAqnpFVQdU9W4A/7FYjaQlgxcuIJFl1JBQxeDFiy63iNx26NQhXJ27im23b8PxruPYdvs2XJ27ikOnDnndtKqTs6uiqi86fWJVfdzp91J2Z2ZmUj3OTNcSCZy9ft3lFpHbzlw6g70P7MXTn3kaNYEajHaN4pk3nsHLP33Z66ZVnRWP9URkCwCo6nuWY78DIGY9RqVxd10d1gcCtkF0fSCAu+rrPWgVuentJ99Oe1wTqMG+B/dh34P7vGlQFXOykb4fQCjjWCOAvtU3h/KJbNqEQJY5roAIIk1NLreIqHo5CaDtqvr31gOq+h0A4eI0iXJpqK3Fay0taKipwfpA8p9vfSCAhpoavNbSwgUkIhc5ebeJiDSo6nTm8WI0iPLbEQzi/PbtGLx4EWevX8dd9fWINDUxeBK5zMk7bgjJ4fp/Mg+IyAsADherUZTfzbW16Nq82etmEFW1FQdQVe0RkXERmUJyQ33I+PxQsRtHRORnjsZ8qtomIg/BCJ6qerS4zSIi8j9HAdTYynSnuU9URH4TyUDK++CJqGqseBXeuC9+GECv5XAzAMeb7omIypGTbUy9SM53plbduY2JiKqRkwC6QVWv2BznNiYiqipOAuhREXkMQCqjhYgMIpnmLi8R6RCRsIjsyXK+z/jc7aBtRESucRJAdwP4cwDNIvJ9EbmE5BzoH+b7RhFpBQBVHQEQNx9n6BaRSSS3RhER+daKA6iRsq4dycz0w0jmCW0vcAU+gmQSZiAZIO3mTXerarMRZImIfMvJKvxvisgtxt7PwwC2icifisgtBXx7EMkSIKaNNteE8gzxu0VkTETGLjL3JRF5yMkQ/kUs1UM6DOAJAA+jSNuYVHW/0fvcKCLLeqhGwuZ2VW1vYuYhIvKQk430zar6nlG+IwygUVV/ZdzamU8cS8E3iGRdpRRj4eiSqg4b5zLT5hER+YaTHqg5BA8DeNcy91nINqZBLAXFEIyVexExK6GNYWk1v9l4TETkS04C6IiInAAwAOCbAGDcF5832Klq1Lg+DCBuPgZw1HL+cRHpADBpOU9E5DtOsjE9KSK/i2QAtCYRKSgjvaoO2Bxry3WeiMiPnGZj+k7GY2ZjIipQvrrtrOtePvKWNTbuOloREblFRL7vvFlElWn25AHMvvUcNEtpalXF7FvPYfbkAZdbRk7kDKBGurqHjTuOfjvfkxmB86sAxgE8VaQ2ElUEVYXOTWP+9JBtEDWD5/zpIejcdNYgS/6RdwhvzHl2AHhRRO5EcpU8huSWpCkkN8OHALQanwcAPKyq75as1URlSESw7v4vAADmTw8BANbd/wWISFrwXLO1M3Wc/K2gOVBjX+awce/6QwDuQ3KbEZAMpJeQ7HGOZMnURETIHkSdBE/OpXpvRYtIxrYibi0iWoXMIGoG0pUEz9mTB6Bz01mvN3u0srYB67Z1FfcHoBQn+0CJaJWsQdS0kp4n51L9gQGUyANmkLPKtTpvZQbfNVs7lwVRzqW6y9E+UPLe9MICBi9cwJmZGdxdV4fIpk1oqOU/ZzmwC3LmY6Cwnmgx51LJOb7jytCb8Th2TUwgoYpriQTWBwL40uQkXmtpwY5gMP8TkGey9RCzrc7nUoy5VFodDuHLzPTCAnZNTGB6cRHXEgkAwLVEAtOLi9g1MYGrCwset5CyyTW8zjUsz2U1c6m0eo4CqIhsMcobm49/s8CEyrRKgxcuIJHljZVQxSCTTPuWiEDWNmTtIVqDqKxtKHhByelcKq0e68KXmTMzM6meZ6ZriQTOXr/ucovcceXGFXzqG5/ClRvlvc143baunD1EM4gWsvUos0d78+feXHEPllaHdeHLzN11dVgfsP9nWx8I4K76epdb5I5Xf/Yq3vnlO/jeme953ZRVK2RuM59cc6kMou5hXfgyE9m0CYEsb7CACCIVWubk4KmDyc9vH/S4Jd4rxVwqOeNkFX5VdeFpdRpqa/FaS8uyVfiACF5racHNFbKV6bunv4sfvvfD1OM33n8DAPD6+6/j83/3+dTxB7c8iMe2rjhhWFkrdC4VQMFzqeSMk3fbbiQzyDcbKevakUwu8lAxG0bZ7QgGcX77dgxevIiz16/jrvp6RJqaKiZ4AsD84jxeGHsBC4n0XQWzi7N4/q3nAQC1gVrs+MQOL5rn2JUbV/DpA5/G8a7juPWmWx0/z7ptXTnvdTeDKINnaTnJSH8FQLtRxiMEYD8TKrvv5tpadG3e7HUzSiZybwQtH2vBo99+FB9Of4iZhZnUubraOmxu2IxXPvsK7mm6x8NWrpx1Lvf3Wn5vVc9VjLlUWp3V7AM9gWSRuBNGHlBuY6KiuqfpHox3j2NucS7t+NziHKLd0bILngDncivNinugIrIbRjE562Ek50RritEoItOx94+hfk09ZhZmUkPWuto6HPvgGB755CNeNy8vzuVWNieTZl9FMvfnMJZKHBOVxKFTh3B17ira/1k7/mrXX+GPXvsjjJ0fw6FTh8oigFbqXC4lORnCi6p+TVXfVdUr1o+it46q3plLZ7D3gb0Y7RrF/Xfcj9GuUex9YC/OTJ3xumkFidwbwaknTyHUGEJdbV3aubraOoQaQzj15Ck8fu/jHrWQVkNWukfMqHn0E1X9m9I0qXDt7e06Npa3HD2R5+I34rht/21Y1MXUsRqpwdSeqVWtxlNpiMi4qrbnu85JD3QQwHdEZEpETlg/HDwXUVUw53JrA7WokRrUBmpRv6Yexz445nXTaBWcBNABJDfNfxXA4YwPIrJhzuVuu30bjncdx7bbt+Hq3FUcOnXI66bRKjhZRGpW1Q1FbwlRBTPncp/+zNOoCdRgtGsUz7zxDF7+6cteN41Wwckc6A8A/K6qTpemSYXjHCgRlUKhc6BOeqBHALwnIocBTFpPqOpfOng+IqKy5CSARgC8i2Rt+PssxxUAAygRVQ0n98Ln7dYSUfUqVsKUcsCaSERUVJWU/DqfvD1QEXkWQL+qvmc8/sNs16rqS8VrGhGVI2vClNVmnPK7QobwnTAWjozHT2a5TgEwgBJVmWpOmJI3gKrqXRmPOQdKRCnVnDClaHOgzAdKVJ2qOWGKk7LGy9LgiMg2AENFaRERlZ1KTH5dCCc90I2ZB1T1JJK1kYioSlVjwpSCA6iI/MAoIneriHw/4+MMmFyZyFNXblzBp77xKVy54U1q3mpMmLKSjfRDSJbu2IlkNnqrS2BZYyJPFbNgnROZCVOO/4fj+Mqxr2RNmJKrqmi5cJJM5LCq+mI2mMlEiJY8/K2HcSR2BA+HHsb3f//7nrZl9uQB6Nx01tLKqorZt56DrG3Aum1dHrQwt5IlE/FL8CSqdn7df6mq0LlpzJ9OritnBlEzeM6fHsKarZ1l3RN1UpVzC4A+AEEAaXlBVfU+m2/J/P4OAHEAraq6P8d1e3KdJ6p2ft1/KSJYd/8Xkm3MCKKZwTNbD7VcOMnGZM5/Dq70G0WkFQBUdUREQiLSqqpRm+vCSM61MoASZRG5N4KWj7Xg0W8/ig+nP8TMwkzqXF1tHTY3bMYrn33Fky1E2YJoJQVPwFkADa0iI30EydtCASAGIAxgWQAlosKY+y9v239b2nFz/6WX2ZAyg6gZSCsleALO9oGOiUiDw9cLIn2707I9pUavlCv6RAXy8/5LaxA1VUrwBJwFUDMj/Qsi8qfWjyK1KWfvVkS6RWRMRMYuXrxYpJckKl9+3n9pznlazb71HFa6+yfX86/m/Go5CaDWjPRPWD4iBXxvHEsBMghgynqykN6nqg6oaruqtjc1Na207UQVx9x/Odo1ivvvuB+jXaPY+8BenJladte1qzIXjG7+3JtYs7UT86eHihJEZ08eyPk85uvPnjywqtfJxe2M9INYuuUzBGPzvYgEVTUOICQiISSD7IZsi0xEtOTtJ99Oe1wTqMG+B/dh34P7vGkQlgdPc9iebXXeyfP7YauUk0Ukx1Q1KiLtxip73BIcjwJoU9VhIDlMR7KHSkRlJtdWpWIFUb9slXKyD3QMyeTJyxSyD1RVB2yOtdlcs+w6IvI/EYGsbcgavKzBT9Y2OA5uftgq5aQHarf/MwLgxCrbQkQVYt22rpzDZjP4rTa4eb1Vyskc6Ncyj4nIi3CwsZ6IKle+4FWs4GYGUTN4Au5tlSpKRnpzAagYz0VEtBKl3iqVi5M5ULv9nhuRZ/8mEVGx2S0YmY+B0vdEncyBPmFzLAaAWZpcMr2wgMELF3BmZgZ319UhsmkTGmpd3VBB5LlSb5UqhNv7QGmV3ozHsWtiAglVXEsksD4QwJcmJ/FaSwt2BLnzi6qDG1ulCrGSkh635Kq8KSK/U5wmUTbTCwvYNTGB6cVFXEskAADXEglMLy5i18QEri4s5HkGospgt1XKWtLEDKJrtnauaqtUPgUFUBH5JpK3YV4Wkf+WcW6LiPwAS1mWqEQGL1xAIsvEeEIVg8wNQFVCVbFuW1daz9Ja0sTcQrXu/i+UNON93gAqIruRvP2yDcn73x8XkS7j3J8hOf95GVxEKrkzMzOpnmema4kEzl6/7nKLyCmvC8CVM+s98Nae5cFTB5Of3z6YugfeD/tAuwF0qOp7ACAiTwJ4QUSeMs63GWWNqcRlUbubAAAPg0lEQVTurqvD+kDANoiuDwRwV329B60iJ7wuAFeurPfAv/zRKfxo7VJmzVRJk/f+Hl+8+I+oabwbNR+dKmlJk0ICaMgMnkAqm3wzgG5VfakkrSJbkU2b8KXJSdtzARFEmJ2qbFh7SwyghbMuEM2MvYAXps5iQdM7FLOJBfRPfwRMf4Tac6MlLWlSSAC9bHMszuDpvobaWrzW0rJsFT4ggtdaWnAztzL5ll8LwJUjM4g+AeBTpw4iEv85Ppq/7klJk0LecXarFqXf4k+2dgSDOL99OwYvXsTZ69dxV309Ik1NDJ4+59cCcOXKDKItAF7/f9/GnT9/K+28WyVNCnnXNYtIZqKQRptjBWVjotW7ubYWXZs3e90MWgE/F4ArV2YQHR3vR53U4IYmoBKAiKCutg7HPjiGRz75SEnbUEgAtauMebTYDSGqdH4uAFeOzM30//vqBVzTRWxbux7/vXU3/vMHP8LY+TEcOnXI+wCqqk/lu4aICmMWgJtZmEltw3Grt1RJrHcixWpvwt7P/AX+rO4WJH76HfxDyxP42l278PI/vVzydhQlGxMRFcbPBeCKwY0ib5m3cZ76k3PY99v7UP+v/gRrtnYi8dPv4Kn6IKI9pa8GxJUHIheZBeCe/szTqAnUYLRrFM+88Qxe/mnpe0ulNnvyAHRuOut952bgk7UNju8O8ss98CYGUCIX+bEAXDGoKuLXLuCBN57B6/Mz2PRbT5WkyJtb5UIKxQBKRKsmIjhyyyfwT/MzePXUX+Pfr6krWZE3t8qFFIJzoERUFOY87v/RQFrt91IUeXOrXEg+7IESkSPZ7q46Fn8PT63fhMUf9QE/6sO/Xncrfrf9Sdd6hW5iACUiR3LdXfU/PjgGAKiFYPu6WyoyeAIcwhORQ5F7Izj15CmEGkOoq61LO1cXWIMttetwfPO/xGPrb3OtyJvbGECJyDHz7qq5xbm043OJefz4gb9A++4o1mztTJsTrSQcwluwWBvRyqXdXZVYhACoq1mLt5ruxaMe7M10E6ODgcXaiJwx765qveXj+Fp9EH92PY7or87hW//4LTz66496ssHdLQygSC/WZjKzvu+amMD57duZLo4oC/Puqj233IHA/DX8uP2P8JVjX0m7u8rtDe5ukXKek2hvb9exsbFVP89L58/ji2fPZi2V8dzddzN9HFEB8t1h5PQOJLeJyHghJdy5iAQWayMqFr9scHcLAyiWirXZYbE2IsqGARTJYm2BLH8ZWayNiLJhAMVSsbaGmppUT3R9IICGmhoWayOirBgZDCzWRkQrxehgwWJtRLQSHMITUUlcuXEFn/rGp3DlxhWvm1IyDKBEVBKv/uxVvPPLd/C9M9/zuiklwwBKRCVx8NTB5Oe3D3rcktLhHCgRFUW2BMuvv/86Pv93n08df3DLg3hs62NuN68kGEBtMCsT0crlSrD8/FvPAwBqA7XY8YkdXjSvJBgVMjArE5EzkXsjaPlYCx799qP4cPpDzCzMpM7V1dZhc8NmvPLZV3BP0z0etrK4OAdqYc3KZN4bfy2RwPTiInZNTODqwkKeZyCqblkTLC/OIdodrajgCTCAphm8cAGJLNmpEqoYvHjR5RYRlR8zwXJtoBY1UoPaQC3q19TjmFEnqZK4HkBFpENEwiKyJ8v5sPHR53bbmJWJaPXMBMvbbt+G413Hse32bbg6dzVV9riSuBpARaQVAFR1BEDcfGw5HwbQaZxvzTxfaqvNyjS9sICXzp9H7+QkXjp/HtMc8lMVMhMsj3aN4v477sdo1yj2PrAXZ6bOeN20onM1obLRqzyiqiNGsGxV1f1Zrp1U1eZcz1eshMqm6YUF3DE6mpaZ3tRQU5MzM73d4lNAhItPRGXIrwmVgwAuWR5vtLvIGN73ZDnXLSJjIjJ2schzkk6zMnHxiag6+XIbk6ruF5EhERlT1XjGuQEAA0CyB1rs13aSlamQxScmKSGqPG4H0DiADcbXQQBT1pOWOdIogBiAbgC2Q/xSWmlWJi4+EVUnt4fwgwBCxtchACMAICLmJGEY6QE25mrrHGJJEKLq5GoANXqW5mp73HwM4KjxeQBASES6jeuH3WyfUywJQlSdXJ8DNeYwM4+1GZ/jMOY3y4m5+JRtFZ5Z7YkqE9/ZRcKSIETVh+/uImJJEKLqwnvhiYgcYgAlInKIAZSIyCEGUCIihxhAiYgc4ip8BWNtJ6LS4rupQrG2E1HpcQhfgZhej8gdDKAViLWdiNzBAFqBmF6PyB0MoBWI6fWI3FE1AbSaCr4xvR6RO6piFb7aVqSZXo/IHRX/TrKuSJvM+cFdExM5K22WM6bXIyq9in83VXPBN6bXIyqtip8D5Yo0EZVKxQdQrkgTUalUfADlijQRlUrFB1BzRbqhpibVE10fCKChpoYr0kS0KlURPbgiTUSlUDURhCvSRFRsFT+EJyIqFQZQIiKHGECJiByqmjnQlWApDCIqBKNChmpLPEJEznEIb8FSGES0EuyBWuRLPHLwo4+wLhDg0J6IADCApsmXeORPzp7F2kCAQ3siAsAhfJpciUcAYB7g0J6IUhhALXIlHsmGVS6JqhcDqEW2xCNrcwRV5hQlql6cA81gl3hkZnERT8VitvOjzClKVL0YQG1kJh6ZXljAn7/7ru21zClKVL04hC8Ac4oSkR2+8wvEnKJElInv/hVgTlEisuIQnojIIQZQIiKHXB/Ci0gHgDiAVlXdb3O+2/iyWVV7XW0cEdEKuNoDFZFWAFDVEQBx87HlfBjAiKoOAAgZj4mIfMntIXwEyd4nAMQAZAbIkOVYzHhMRORLbg/hgwAuWR5vtJ40ep6mVgCDbjSKiMgJXy4iGUP7qKpGbc51i8iYiIxdZBIPIvKQ2wE0DmCD8XUQwFSW68LZFpBUdUBV21W1vYm3UBKRh9wOoINYmtcMARgBABFJZSQWkW5zdZ6LSETkZ64GUHNIbgTGuGWIftRyvE9EJkXkspttIyJaKdf3gWYsFJnH2ozPIwAa3W4TEZETvlxEIiIqBwygREQOiWYp41sOROQigPeL9HS3AfhlkZ6rFPzcPrbNGbbNGTfa9s9VNe82n7IOoMUkImOq2u51O7Lxc/vYNmfYNmf81DYO4YmIHGIAJSJyiAF0ybLtVT7j5/axbc6wbc74pm2cAyUicog9UCIih6o6gNokdG4VkQ4ja77nrO0z2qbGba6TItLvl7YZjztEJGypKOAZm7btMdrnedv8yshy1i0ifRnHW7N9D1VxADXuux/KOPxlVR1GMhu+p/9xbNq3QVVFVZsBdALos//O0stsm/G7ihm34sa8/N3ZtC0MAMa/a7OIeJak2y5IWf7w7PGwXbaVILK8R1yX5fdmG/DdVrUB1Hyzm4+NXucJ49x+u1ykbspsn/HY1K6qseXf5Y7MthnM/8ghL393Nm3baXk8ieVVEFxhF6TylbhxkW0liCz/zq7K8nvzTemfqg2gNu4DsNEYKnvWG8jH+M9y2Ot2WBkBM2Zk0LqU73qXTSE9B22zR+2wC1L5Sty4wsixa65stwIY86IdWdj93nxT+ocBNN2UJeWeL+ZBbexU1Xj+y9xj5HONA3gWwIteDpNtDGMpaDYjexLvksoSpHKWuHFbrkoQXrH7vfkp4DOALpnC0nAljmSP1I/8OKnfDeBZIxH2bgC++eNjTHUMGsEhDu+HpL4LUhZZK0F4ze735offJQPokmEsDQWCMOZD/cRnPTtbxmKNb3rIxpus3XiTBY32eckapAotcVNyZVAJwi64ex7wqzaAGkP0dnOobvRU4sbjjV6/0TLbZ+FpDwqw/d3tB9BtbhWyS5rtYduiAC4Zj73e+pUZpGxL3HjQLttKEDn+D7rKLrj7JeDzTiQiF1i2BF1CstfZqaojxt7UGJK7F3xzi6Jf2P3ejFPLfpeetI8BlIjImaodwhMRrRYDKBGRQwygREQOMYASETnEAEpUYqXIAsXMUv7AAEoAkpv0RWTI2Atops3rM27TzPe9k8XIilPE59mT8XP0W38O4+dSEVmWaUhExq25ECzXmmkELxvXFBTAjLSDWfMDGG09UsDrZe7FvGTXfnIXAyiZG6Ynkbz7qk1VBcn9doXeNtqD4mxSX/XzGEElgmTOAPPn2GBzaQxAR4F3d0VVtdn4aDSes0fy5GQ1fq8hBzdlZL5eL5I5BlKB3XjOIHui3uI+0Cpn9MwuA+iphI3cIqIAmnOl+zN6uanAqaqdlnPjAAYtd7n0IXnLYFvGc7QCGAfQmC25i4hMIvl7zbrJ2wiKO1V1Z57X6wDwohFQrW0YMnLEkgfYA6U+JJMhl33wtCi05/wsCu+FrogR3EJFvENmWS/auE11g0/vXa8KDKDUjmQilZyM7N/mnNyk9U2bOSdoziOKyBHLHF7eIFWk59kPYMiY9+zINYdrBKARJIfIBbPcXrg/R2rBMGzubTfmmseNn+cI8mT9EhFzmN6HpdsYrUaQTBpNHmAApRAKywJ0Ccn50UYk38z5FjC+jOQ9yuaQ0+ni0Iqex8jO04nkzzUE4LLkTpDdi2QilFyBudUIeJeNKYJ+AH15MgHdB/vEL0eQzGnZaAzb7QJ86vWQnF7pQ/b7vWPwZ4rDqsAASjFkJPI1A4XxYWY1GjZ7W8ZwPyi5S1ActvTOBmEfKAqx4ucx2mpdROrLtthSYC80agS8RiQXulDAlIeZZDrF6LmGVLXHcvgIlku9XsbPYLdoZc24Ty5jAKUxZCRAtrxxUwHAGEr2GcPP8QKedzLjsdM3+aqex1it3g/74a/J7IUWEpwHgNTiTz6ZPfsQgBUn/1XVEWNRqdtmO5Nvcq9WIwZQ6kWyMFfW7TBGYHkXwAlVbctcIS4Dy3qDVpZeaKHTDGaPMNewP47lJTpWW7/HmoDZ1Az/1aGqGgygVc4YHncC6DcWbILAsnrgG5AMQr6uFyXJgoCXjQWvkGUBphvJFfdceo3rCumFRpFceMu1D/QSMoKlOYdpDsWN33XP8m9NZ/wcZv7LzIKCQfggyXa1YgClVM10JBc+xo2FkiEk36wjxp7KYQCTxt5GX9aLMgJbp/ExDmN/K5L7LHMOnS290EJ7iLsBmCV27RyB/eJOG5JZ3i8DeBH2WehbrXciIfmzmIt4mT3pduM8eYAb6YlKwHKDQtaN9uXwGpQbe6BEJWAEtAEkt2GVyuMAhhk8vcMeKFGJGItMR2A/9C7G808iOT3BOVCPsAdKVCJGYOtE/psOVsxYVOph8PQWe6BERA6xB0pE5BADKBGRQwygREQOMYASETnEAEpE5BADKBGRQ/8fsMimSkkYAmAAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 360x360 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig = plt.figure(figsize=(5,5))\n",
    "plt.scatter(np.array(l1_snrs)-INPUT_SNR, l1_times, label='vector-GTF+L1', s=50, c='c')\n",
    "plt.scatter(np.array(scad_snrs)-INPUT_SNR, np.array(l1_times)+np.array(scad_times), label='vector-GTF+SCAD', marker='x', s=100, c='#f29f3a')\n",
    "plt.scatter(np.array(mcp_snrs)-INPUT_SNR, np.array(l1_times)+np.array(mcp_times), label='vector-GTF+MCP', marker='*',s=100,c='g')\n",
    "plt.xlabel('Gain in SNR (dB)',fontsize=16)\n",
    "plt.ylabel('Runtime (sec)',fontsize=16)\n",
    "plt.legend(fontsize='large')\n",
    "fig.savefig('2d-grid-runtime.pdf',format = 'pdf',bbox_inches='tight')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Minnesota runtime"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    input_snr (dB)  sigma_sq penalty_f  gamma_mult  rho_mult  penalty_param  \\\n",
      "0              -10   44.1938        L1      0.7725    0.1636            0.0   \n",
      "1              -10   44.1938        L1      0.4085    0.1567            0.0   \n",
      "2              -10   44.1938   SCAD-L1     63.1302   10.8704            3.7   \n",
      "3              -10   44.1938   SCAD-L1     63.4909   19.6460            3.7   \n",
      "4              -10   44.1938    MCP-L1    146.0735   19.7557            1.4   \n",
      "5              -10   44.1938    MCP-L1     93.4036   19.5608            1.4   \n",
      "6               -5   13.9753        L1      1.3504    4.1660            0.0   \n",
      "7               -5   13.9753        L1      0.6381    1.3307            0.0   \n",
      "8               -5   13.9753   SCAD-L1      2.5162    0.7747            3.7   \n",
      "9               -5   13.9753   SCAD-L1    142.2719   19.5464            3.7   \n",
      "10              -5   13.9753    MCP-L1    130.2446    2.5336            1.4   \n",
      "11              -5   13.9753    MCP-L1    124.4601   19.3936            1.4   \n",
      "12               0    4.4194        L1      2.2541   17.9244            0.0   \n",
      "13               0    4.4194        L1      1.0972   18.3523            0.0   \n",
      "14               0    4.4194   SCAD-L1      2.5999    0.4941            3.7   \n",
      "15               0    4.4194   SCAD-L1    116.2501   20.0837            3.7   \n",
      "16               0    4.4194    MCP-L1      3.1983    0.2482            1.4   \n",
      "17               0    4.4194    MCP-L1    139.1195   16.1958            1.4   \n",
      "18               5    1.3975        L1      3.9740    7.9703            0.0   \n",
      "19               5    1.3975        L1      1.7572    0.4986            0.0   \n",
      "20               5    1.3975   SCAD-L1      7.8043    1.3317            3.7   \n",
      "21               5    1.3975   SCAD-L1    141.8960   17.7472            3.7   \n",
      "22               5    1.3975    MCP-L1     13.2344    0.7882            1.4   \n",
      "23               5    1.3975    MCP-L1      9.5121    1.8951            1.4   \n",
      "24              10    0.4419        L1      6.2723   13.2339            0.0   \n",
      "25              10    0.4419        L1      2.8601    9.8669            0.0   \n",
      "26              10    0.4419   SCAD-L1    141.3012    2.0873            3.7   \n",
      "27              10    0.4419   SCAD-L1      5.6829    5.4911            3.7   \n",
      "28              10    0.4419    MCP-L1     38.3344    0.8567            1.4   \n",
      "29              10    0.4419    MCP-L1      9.2113    2.2030            1.4   \n",
      "30              15    0.1398        L1     10.6506    8.7773            0.0   \n",
      "31              15    0.1398        L1      4.5444   18.0948            0.0   \n",
      "32              15    0.1398   SCAD-L1     20.2751    2.6592            3.7   \n",
      "33              15    0.1398   SCAD-L1     10.0412    6.6922            3.7   \n",
      "34              15    0.1398    MCP-L1     19.9688    1.0773            1.4   \n",
      "35              15    0.1398    MCP-L1      7.8069    1.4974            1.4   \n",
      "36              20    0.0442        L1     18.9923    8.9085            0.0   \n",
      "37              20    0.0442        L1      7.3462    0.5452            0.0   \n",
      "38              20    0.0442   SCAD-L1    133.0239    4.2697            3.7   \n",
      "39              20    0.0442   SCAD-L1     26.8969    6.7859            3.7   \n",
      "40              20    0.0442    MCP-L1    146.5169    1.5199            1.4   \n",
      "41              20    0.0442    MCP-L1     35.4461    2.7203            1.4   \n",
      "42              25    0.0140        L1     36.2243   16.6506            0.0   \n",
      "43              25    0.0140        L1     14.1000   13.7031            0.0   \n",
      "44              25    0.0140   SCAD-L1    139.9638    7.8425            3.7   \n",
      "45              25    0.0140   SCAD-L1     47.9893    9.7751            3.7   \n",
      "46              25    0.0140    MCP-L1    147.4437    2.2780            1.4   \n",
      "47              25    0.0140    MCP-L1     98.3396    4.6010            1.4   \n",
      "48              30    0.0044        L1     60.4261   19.2467            0.0   \n",
      "49              30    0.0044        L1     25.3854   19.0669            0.0   \n",
      "50              30    0.0044   SCAD-L1    132.1196   12.7539            3.7   \n",
      "51              30    0.0044   SCAD-L1    141.2749   16.8143            3.7   \n",
      "52              30    0.0044    MCP-L1     96.5133    4.8295            1.4   \n",
      "53              30    0.0044    MCP-L1    145.5949    9.1012            1.4   \n",
      "\n",
      "    avg output_snr (dB)        0        1        2        3        4        5  \\\n",
      "0               10.7008  10.8650  10.4027  10.4276  10.0994   9.3455  11.2891   \n",
      "1                8.9006   8.5979   8.1463   8.8126   8.3962   7.5783   9.8721   \n",
      "2               10.6510  11.1753  10.6343  10.4650  10.3519   9.2945  11.2464   \n",
      "3                8.5527   8.3265   7.9038   8.5624   8.4977   7.2399   9.4344   \n",
      "4               10.6934  11.0866  10.5421  10.3535  10.1414   9.3022  11.2920   \n",
      "5                8.7515   8.5095   7.9915   8.6623   8.3514   7.5471   9.4495   \n",
      "6               13.8589  12.7627  13.1180  14.0752  14.2906  13.6382  14.8464   \n",
      "7               11.5855  10.8032  11.0432  11.7642  11.8061  11.4952  12.4575   \n",
      "8               15.0252  14.1109  15.0644  14.2583  15.4779  14.7800  15.5852   \n",
      "9               11.2520  10.5680  10.4820  11.6576  11.9091  10.8956  11.5289   \n",
      "10              14.4758  13.9652  13.0991  13.4198  15.8513  14.1935  14.7184   \n",
      "11              11.5180  10.7471  11.1864  11.7356  12.0007  11.2767  12.0760   \n",
      "12              16.6669  16.8474  17.3664  16.4441  16.8156  16.6465  15.8054   \n",
      "13              14.1629  14.4290  15.0236  13.8498  14.2058  14.2324  13.3097   \n",
      "14              17.9338  18.0289  18.3659  18.1752  18.0514  17.6783  17.3114   \n",
      "15              14.1378  14.6763  15.4008  13.7299  13.8414  13.8127  13.2286   \n",
      "16              18.3211  18.3783  19.0201  18.2270  17.7523  19.2688  17.3555   \n",
      "17              14.0021  14.2557  15.0730  13.5305  13.6834  13.9785  13.2115   \n",
      "18              20.4455  20.2792  21.9008  19.9322  20.2176  20.7260  20.3626   \n",
      "19              17.6019  17.3292  19.2596  17.1332  17.1147  17.9370  17.3772   \n",
      "20              23.3067  23.9769  23.0736  23.7184  23.4573  24.8348  23.3608   \n",
      "21              17.9980  17.9587  19.9677  17.8139  17.5122  18.8434  17.6421   \n",
      "22              23.2318  24.4798  23.0012  23.3124  23.3007  25.4584  23.1263   \n",
      "23              17.5001  17.1457  19.1241  15.9306  16.8126  18.6359  17.4525   \n",
      "24              24.8654  25.0759  25.7643  24.7262  24.2799  25.0607  24.7445   \n",
      "25              21.7610  21.9395  22.9600  21.3233  21.0903  22.2484  21.4725   \n",
      "26              29.8238  29.9879  29.7902  30.2575  30.3172  29.8820  30.1367   \n",
      "27              23.1640  23.3194  25.5453  22.1933  22.0261  24.0651  22.5582   \n",
      "28              29.4451  29.3593  29.6534  30.1034  29.2010  29.3713  29.3556   \n",
      "29              23.1255  23.3206  25.4995  22.1431  21.5680  24.1490  23.6626   \n",
      "30              29.2096  29.3316  29.3400  29.4732  28.6106  29.3176  29.6270   \n",
      "31              25.9752  26.1764  26.1993  25.9494  25.4717  26.0906  26.2838   \n",
      "32              38.2706  39.2370  39.8538  38.9752  36.9283  35.9711  39.1723   \n",
      "33              29.6625  31.1278  32.1340  29.0893  27.7982  27.8205  29.2116   \n",
      "34              37.5236  38.4008  38.6399  38.5204  35.5103  35.8147  38.2301   \n",
      "35              30.0238  31.6634  32.4174  30.4481  28.8477  27.8109  29.1072   \n",
      "36              34.0417  33.8525  33.7978  34.2493  33.4885  34.2598  33.9488   \n",
      "37              30.3081  30.3577  30.1591  30.2311  29.7585  30.9084  30.2870   \n",
      "38              47.7988  47.7988  50.4533  49.9640  48.6376  48.1023  43.0017   \n",
      "39              38.0247  37.2396  37.4269  39.8683  39.9998  37.1208  38.6479   \n",
      "40              47.2294  47.5423  49.9859  48.2860  48.3511  47.0399  42.9685   \n",
      "41              37.9404  39.4729  37.1843  39.5703  39.7029  36.6599  38.3432   \n",
      "42              39.2230  38.2002  39.5485  38.9143  38.9734  40.0243  38.6878   \n",
      "43              35.9132  35.1251  36.1296  35.9057  35.4016  36.4792  35.8914   \n",
      "44              54.1011  60.5911  51.9721  58.7920  55.7790  48.4265  55.1848   \n",
      "45              51.6040  49.7889  50.5452  52.4196  54.7542  47.8436  52.3819   \n",
      "46              50.9520  48.3622  50.3307  52.4056  54.1130  47.5429  50.8442   \n",
      "47              51.7811  50.1060  50.4509  52.3896  54.7926  47.9284  52.8527   \n",
      "48              44.0462  44.4073  43.4485  44.2419  44.2523  44.2726  44.8821   \n",
      "49              40.6292  40.9299  39.9856  40.3030  41.0067  40.8022  41.2378   \n",
      "50              57.1988  59.6659  56.5018  55.7028  58.9472  56.1746  59.1584   \n",
      "51              57.3008  61.4285  56.7578  55.7960  59.3595  56.4705  59.0485   \n",
      "52              54.5998  58.1200  52.3865  54.7270  55.2140  53.5633  56.5649   \n",
      "53              53.4034  53.3520  55.2539  50.4290  52.6663  53.0776  57.3214   \n",
      "\n",
      "          6        7        8        9  \n",
      "0   11.1417  10.7134  11.6011  11.6717  \n",
      "1    9.4529   9.3520   9.4033  10.0560  \n",
      "2   10.7998  10.5290  10.7167  11.7493  \n",
      "3    8.9692   8.6998   8.8648   9.5232  \n",
      "4   11.0295  10.6761  11.3531  11.6802  \n",
      "5    9.0599   9.1977   9.1703  10.1705  \n",
      "6   14.0295  14.7783  13.3371  14.2069  \n",
      "7   11.4118  12.6282  11.1100  11.6855  \n",
      "8   15.3997  15.5731  14.8654  15.4492  \n",
      "9   11.2026  12.4084  11.3175  10.9236  \n",
      "10  15.0334  15.0668  15.0611  15.1343  \n",
      "11  11.1067  12.4299  11.2763  11.6196  \n",
      "12  17.7593  16.4032  16.1847  16.7166  \n",
      "13  14.9616  13.8818  13.7417  14.2846  \n",
      "14  18.8297  17.8092  17.3355  17.9685  \n",
      "15  15.1563  14.0583  13.6079  14.3597  \n",
      "16  18.9772  18.2259  18.1330  18.2317  \n",
      "17  14.4357  13.4335  14.2073  14.5538  \n",
      "18  20.6699  19.2746  21.0385  20.5409  \n",
      "19  17.9208  16.5525  18.1634  17.7809  \n",
      "20  23.4398  23.1033  22.2613  22.4021  \n",
      "21  18.3670  16.7392  17.9960  17.8693  \n",
      "22  23.6464  23.0783  21.6487  22.3645  \n",
      "23  19.4858  17.4016  17.1429  17.0960  \n",
      "24  23.5922  25.3741  25.1073  25.3363  \n",
      "25  20.6629  22.0419  21.9263  22.4252  \n",
      "26  28.9259  30.9028  30.0526  28.5003  \n",
      "27  23.5178  23.8124  22.9619  22.6948  \n",
      "28  28.9130  30.8494  29.9602  28.2079  \n",
      "29  22.8100  23.5751  22.7950  22.9135  \n",
      "30  28.7245  29.0589  28.9859  29.7736  \n",
      "31  25.6598  25.6950  25.8800  26.4451  \n",
      "32  38.8006  37.8917  37.8883  39.7348  \n",
      "33  29.3423  29.5104  30.9367  32.2778  \n",
      "34  38.0409  37.7287  37.6191  38.1130  \n",
      "35  30.4203  28.7478  30.1970  33.9426  \n",
      "36  34.0405  34.0995  34.6200  34.1585  \n",
      "37  30.0905  30.2737  30.5992  30.5158  \n",
      "38  52.2601  47.2362  51.4110  47.2164  \n",
      "39  37.3875  39.6555  39.9345  35.5235  \n",
      "40  51.4203  46.4044  49.9321  46.6448  \n",
      "41  37.1064  39.1558  39.5292  35.2965  \n",
      "42  39.6137  39.6002  39.5095  39.4726  \n",
      "43  36.4606  36.1938  36.0932  35.6538  \n",
      "44  53.6492  55.5476  57.5375  58.4045  \n",
      "45  52.7754  53.4329  53.4097  53.8268  \n",
      "46  52.6397  52.5386  52.9351  53.1023  \n",
      "47  53.0637  53.4449  53.8934  54.3655  \n",
      "48  43.2441  43.6079  44.5529  43.8362  \n",
      "49  40.0940  40.5344  40.9750  40.6056  \n",
      "50  54.4881  56.7379  58.8605  59.3035  \n",
      "51  54.6922  56.5385  57.5878  59.4266  \n",
      "52  53.1338  55.0453  55.1568  54.7671  \n",
      "53  50.4216  55.2698  55.2145  57.5337  \n"
     ]
    }
   ],
   "source": [
    "data = pd.read_csv('../../public_code/datasets/minnesota/minnesota-simulation-19-01-09-03-56.csv',\n",
    "                       skiprows=3,\n",
    "                      names=['input_snr (dB)', 'sigma_sq', 'penalty_f', 'gamma_mult', 'rho_mult', 'penalty_param',\n",
    "                              'avg output_snr (dB)'] + [str(e) for e in range(10)])\n",
    "print(data)\n",
    "\n",
    "input_snr = np.array(sorted(list(set(data['input_snr (dB)'].values))))\n",
    "\n",
    "L1 = data.loc[data['penalty_f'] == 'L1', ['gamma_mult', 'rho_mult']].values\n",
    "SCAD = data.loc[data['penalty_f'] == 'SCAD-L1', ['gamma_mult', 'rho_mult']].values\n",
    "MCP = data.loc[data['penalty_f'] == 'MCP-L1', ['gamma_mult', 'rho_mult']].values\n",
    "    \n",
    "L1_vec = L1[range(0,len(L1),2)]\n",
    "SCAD_vec = SCAD[range(0,len(SCAD),2)]\n",
    "MCP_vec = MCP[range(0,len(MCP),2)]\n",
    "    \n",
    "L1 = L1[range(1,len(L1),2)]\n",
    "SCAD = SCAD[range(1,len(SCAD),2)]\n",
    "MCP = MCP[range(1,len(MCP),2)]\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "n = 2642 m = 3304\n",
      "sparsity 83\n"
     ]
    }
   ],
   "source": [
    "name = 'minnesota'\n",
    "G = pg.graphs.Minnesota(connect=True) \n",
    "Gnx = nx.from_scipy_sparse_matrix(G.A.astype(float), edge_attribute='weight') \n",
    "y_true = np.load('../../public_code/datasets/minnesota/beta_0.npy')\n",
    "n = nx.number_of_nodes(Gnx)\n",
    "k = 0\n",
    "Dk = penalty_matrix(Gnx, k)\n",
    "print ('n =', n,  'm =', Gnx.number_of_edges())\n",
    "print ('sparsity', np.count_nonzero(Dk*y_true))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "DTD = Dk.T.dot(Dk).toarray()\n",
    "[S, V] = np.linalg.eigh(DTD)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "20\n",
      "sigma squared: 0.04419379258137775\n",
      "input signal MSE: 116.28172270773021\n",
      "input signal SNR: 20.017826291195547\n",
      "L1 vec MSE: 5.669266964805986\n",
      "L1 vec SNR: 33.13767177991045\n",
      "100%|██████████| 30/30 [05:04<00:00, 27.09s/it, best loss: -46.78086640516967] \n",
      "100%|██████████| 30/30 [06:41<00:00, 24.63s/it, best loss: -43.90870256285072] \n",
      "sigma squared: 0.04419379258137775\n",
      "input signal MSE: 116.08327045539836\n",
      "input signal SNR: 20.02524451180536\n",
      "L1 vec MSE: 5.392083729944801\n",
      "L1 vec SNR: 33.3553745887724\n",
      "100%|██████████| 30/30 [04:22<00:00,  7.13s/it, best loss: -45.97396938544265]\n",
      "100%|██████████| 30/30 [08:13<00:00, 17.20s/it, best loss: -43.904481972752045]\n",
      "sigma squared: 0.04419379258137775\n",
      "input signal MSE: 117.50275968270603\n",
      "input signal SNR: 19.9724601970059\n",
      "L1 vec MSE: 5.740300397970344\n",
      "L1 vec SNR: 33.083594660839914\n",
      "100%|██████████| 30/30 [03:12<00:00,  9.61s/it, best loss: -48.14785197700323] \n",
      "100%|██████████| 30/30 [07:16<00:00, 25.25s/it, best loss: -38.764425594627184]\n",
      "sigma squared: 0.04419379258137775\n",
      "input signal MSE: 117.31259754879481\n",
      "input signal SNR: 19.979494352220676\n",
      "L1 vec MSE: 5.7141678755393155\n",
      "L1 vec SNR: 33.10341091020146\n",
      "100%|██████████| 30/30 [05:00<00:00, 13.51s/it, best loss: -47.70647420821571] \n",
      "100%|██████████| 30/30 [05:47<00:00, 14.57s/it, best loss: -47.03811835646353]\n",
      "sigma squared: 0.04419379258137775\n",
      "input signal MSE: 116.9283814117695\n",
      "input signal SNR: 19.99374148340623\n",
      "L1 vec MSE: 5.523016394010609\n",
      "L1 vec SNR: 33.25117753970848\n",
      "100%|██████████| 30/30 [03:53<00:00,  7.26s/it, best loss: -47.1029750276221]  \n",
      "100%|██████████| 30/30 [05:25<00:00, 27.24s/it, best loss: -44.843583300628104]\n",
      "sigma squared: 0.04419379258137775\n",
      "input signal MSE: 117.18595329107941\n",
      "input signal SNR: 19.984185290176658\n",
      "L1 vec MSE: 5.494529043411899\n",
      "L1 vec SNR: 33.27363612961901\n",
      "100%|██████████| 30/30 [04:03<00:00, 10.26s/it, best loss: -48.88477507439896]\n",
      "100%|██████████| 30/30 [04:09<00:00, 10.04s/it, best loss: -42.351755925393526]\n",
      "sigma squared: 0.04419379258137775\n",
      "input signal MSE: 115.15401223970795\n",
      "input signal SNR: 20.060150118992826\n",
      "L1 vec MSE: 5.399118733007213\n",
      "L1 vec SNR: 33.349712080903615\n",
      "100%|██████████| 30/30 [04:03<00:00,  9.29s/it, best loss: -47.022862855041026]\n",
      "100%|██████████| 30/30 [07:07<00:00, 26.44s/it, best loss: -40.36277401692796] \n",
      "sigma squared: 0.04419379258137775\n",
      "input signal MSE: 116.23983289732817\n",
      "input signal SNR: 20.019391093534416\n",
      "L1 vec MSE: 5.285258470107026\n",
      "L1 vec SNR: 33.44227855430074\n",
      "100%|██████████| 30/30 [02:32<00:00,  7.53s/it, best loss: -46.364763159807026]\n",
      "100%|██████████| 30/30 [06:28<00:00, 17.58s/it, best loss: -40.00047403299391]\n",
      "sigma squared: 0.04419379258137775\n",
      "input signal MSE: 117.63316873834712\n",
      "input signal SNR: 19.967642903498916\n",
      "L1 vec MSE: 5.403949766704615\n",
      "L1 vec SNR: 33.345827830023254\n",
      "100%|██████████| 30/30 [03:00<00:00, 11.21s/it, best loss: -46.764858207535525]\n",
      "100%|██████████| 30/30 [04:45<00:00, 17.55s/it, best loss: -44.188949801863316]\n",
      "sigma squared: 0.04419379258137775\n",
      "input signal MSE: 116.31515550999556\n",
      "input signal SNR: 20.016577806542646\n",
      "L1 vec MSE: 5.243779326559557\n",
      "L1 vec SNR: 33.47649679284584\n",
      "100%|██████████| 30/30 [02:47<00:00, 18.63s/it, best loss: -47.54952700365628] \n",
      "100%|██████████| 30/30 [06:06<00:00,  9.80s/it, best loss: -45.19294690220932] \n"
     ]
    }
   ],
   "source": [
    "num_trial = 10\n",
    "max_evals = 30\n",
    "d = 20\n",
    "pspace = (\n",
    "    hp.loguniform('gamma_mult', -3, 5),\n",
    "    hp.loguniform('rho_mult', -3, 3)\n",
    ")\n",
    "\n",
    "idx = 6\n",
    "\n",
    "INPUT_SNR = input_snr[idx]\n",
    "print(INPUT_SNR)\n",
    "\n",
    "Ys = []\n",
    "scad_opts = []\n",
    "mcp_opts = []\n",
    "\n",
    "y_true_norm_sq = np.linalg.norm(y_true, 2) ** 2\n",
    "sigma_sq = y_true_norm_sq/(10**(INPUT_SNR/10.0))/n\n",
    "\n",
    "for trial in range(num_trial):\n",
    "    Y = np.tile(y_true.reshape(-1, 1), (1, d)) + np.random.normal(scale=np.sqrt(sigma_sq), size=(n, d))\n",
    "    Ys.append(Y)\n",
    "    ses = np.linalg.norm(Y - np.tile(y_true.reshape(-1, 1), (1, d)), 2, axis=0) ** 2\n",
    "\n",
    "    mse = np.mean(ses)\n",
    "    mean_snr = 10*np.log10(y_true_norm_sq/mse)\n",
    "    \n",
    "    print ('sigma squared:', sigma_sq)\n",
    "    print ('input signal MSE:', mse)\n",
    "    print ('input signal SNR:', mean_snr)\n",
    "\n",
    "# opt_param, B_hat, snr_avg, penalty_param, output_snr = autotune_denoising(Gnx, Y, y_true.reshape(-1,1), 0,\n",
    "#                                                                                           'L1', sigma_sq,\n",
    "#                                                                                           max_evals, pspace,\n",
    "#                                                                                           B_init=None, vec=True)\n",
    "\n",
    "    res = admm_denoising_snr(L1_vec[idx], Y, sigma_sq, y_true, Dk, 'L1', eig=(S, V), B_init=None, vec=True)\n",
    "\n",
    "    B_hat = res['B']\n",
    "    l1_snr = -res['loss']\n",
    "    mse = res['mse']\n",
    "    \n",
    "    print ('L1 vec MSE:', mse)\n",
    "    print ('L1 vec SNR:', l1_snr)\n",
    "\n",
    "\n",
    "    l1_B = B_hat.copy()\n",
    "\n",
    "    opt_param, B_hat, snr_avg, penalty_param, output_snr = autotune_denoising(Y, y_true.reshape(-1,1), Dk,\n",
    "                                                                                          'SCAD', sigma_sq,\n",
    "                                                                                          max_evals, pspace,eig=(S, V),\n",
    "                                                                                          B_init=l1_B, vec=True)\n",
    "    scad_opts.append(opt_param.copy())\n",
    "     \n",
    "    opt_param, B_hat, snr_avg, penalty_param, output_snr = autotune_denoising(Y, y_true.reshape(-1,1), Dk,\n",
    "                                                                                          'MCP', sigma_sq,\n",
    "                                                                                          max_evals, pspace,eig=(S, V),\n",
    "                                                                                          B_init=l1_B, vec=True)\n",
    "    mcp_opts.append(opt_param.copy())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[{'gamma_mult': 140.4774489307071, 'rho_mult': 1.6213601463049063},\n",
       " {'gamma_mult': 136.43107850030668, 'rho_mult': 1.4946109050922964},\n",
       " {'gamma_mult': 146.79815940940836, 'rho_mult': 2.593740505092685},\n",
       " {'gamma_mult': 117.68983206081035, 'rho_mult': 2.898564835803327},\n",
       " {'gamma_mult': 145.2083726635658, 'rho_mult': 3.4406907958173925},\n",
       " {'gamma_mult': 80.01140335584057, 'rho_mult': 3.105911538008845},\n",
       " {'gamma_mult': 76.8992951175442, 'rho_mult': 3.9444993814374727},\n",
       " {'gamma_mult': 146.66575230614012, 'rho_mult': 1.4215590632021096},\n",
       " {'gamma_mult': 141.1612320219524, 'rho_mult': 2.097427157831775},\n",
       " {'gamma_mult': 135.0284499818782, 'rho_mult': 1.8943562564983663}]"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "[{'gamma_mult': 141.27143414687043, 'rho_mult': 0.6900637878333573},\n",
       " {'gamma_mult': 80.99869203700207, 'rho_mult': 1.2945946792672984},\n",
       " {'gamma_mult': 45.26177790637451, 'rho_mult': 0.3139051491806936},\n",
       " {'gamma_mult': 139.59107806815706, 'rho_mult': 1.2687094885406511},\n",
       " {'gamma_mult': 43.63817916509491, 'rho_mult': 1.2947653972972517},\n",
       " {'gamma_mult': 62.768591754739994, 'rho_mult': 0.8739865030579489},\n",
       " {'gamma_mult': 28.999084817746855, 'rho_mult': 1.0993674505816262},\n",
       " {'gamma_mult': 37.20048104582307, 'rho_mult': 0.8922786314602457},\n",
       " {'gamma_mult': 62.59088556360559, 'rho_mult': 1.586730818081994},\n",
       " {'gamma_mult': 146.13964976480858, 'rho_mult': 1.7445244625396426}]"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "sigma squared: 0.04419379258137775\n",
      "input signal MSE: 116.31515550999556\n",
      "input signal SNR: 20.016577806542646\n",
      "L1 0.8393417668433006 7.477276129923544 0\n",
      "L1\n",
      "output signal MSE: 5.669266964805986\n",
      "output signal SNR: 33.13767177991045\n",
      "iteration: (52,)\n",
      "SCAD 6.208231240404755 10.065778712237345 3.7\n",
      "SCAD\n",
      "output signal MSE: 0.24502333447565028\n",
      "output signal SNR: 46.78086640516967\n",
      "iteration: (38,)\n",
      "MCP 6.243320458360557 4.3082893641537785 1.4\n",
      "MCP\n",
      "output signal MSE: 0.47470502615301813\n",
      "output signal SNR: 43.90870256285072\n",
      "iteration: (22,)\n",
      "L1 0.8393417668433006 7.477276129923544 0\n",
      "L1\n",
      "output signal MSE: 5.392083729944801\n",
      "output signal SNR: 33.3553745887724\n",
      "iteration: (52,)\n",
      "SCAD 6.029406784896219 9.011617131943371 3.7\n",
      "SCAD\n",
      "output signal MSE: 0.29505103863225796\n",
      "output signal SNR: 45.97396938544265\n",
      "iteration: (35,)\n",
      "MCP 3.579639395246163 4.634182114781293 1.4\n",
      "MCP\n",
      "output signal MSE: 0.47516658145173807\n",
      "output signal SNR: 43.904481972752045\n",
      "iteration: (32,)\n",
      "L1 0.8393417668433006 7.477276129923544 0\n",
      "L1\n",
      "output signal MSE: 5.740300397970344\n",
      "output signal SNR: 33.083594660839914\n",
      "iteration: (53,)\n",
      "SCAD 6.487567408267419 16.827066366342375 3.7\n",
      "SCAD\n",
      "output signal MSE: 0.17885821352751044\n",
      "output signal SNR: 48.14785197700323\n",
      "iteration: (52,)\n",
      "MCP 2.000289624658701 0.7142857142857143 1.4\n",
      "MCP\n",
      "output signal MSE: 1.5518563829302123\n",
      "output signal SNR: 38.764425594627184\n",
      "iteration: (25,)\n",
      "L1 0.8393417668433006 7.477276129923544 0\n",
      "L1\n",
      "output signal MSE: 5.7141678755393155\n",
      "output signal SNR: 33.10341091020146\n",
      "iteration: (53,)\n",
      "SCAD 5.201160027032633 15.075899559742671 3.7\n",
      "SCAD\n",
      "output signal MSE: 0.19799155451929998\n",
      "output signal SNR: 47.70647420821572\n",
      "iteration: (79,)\n",
      "MCP 6.169059150355041 7.826743879423968 1.4\n",
      "MCP\n",
      "output signal MSE: 0.23093100770037553\n",
      "output signal SNR: 47.03811835646353\n",
      "iteration: (60,)\n",
      "L1 0.8393417668433006 7.477276129923544 0\n",
      "L1\n",
      "output signal MSE: 5.523016394010609\n",
      "output signal SNR: 33.25117753970848\n",
      "iteration: (52,)\n",
      "SCAD 6.417308702573029 22.079974986861874 3.7\n",
      "SCAD\n",
      "output signal MSE: 0.22750795334286772\n",
      "output signal SNR: 47.1029750276221\n",
      "iteration: (115,)\n",
      "MCP 1.9285366386512044 2.497002507145533 1.4\n",
      "MCP\n",
      "output signal MSE: 0.38276811750658546\n",
      "output signal SNR: 44.843583300628104\n",
      "iteration: (501,)\n",
      "L1 0.8393417668433006 7.477276129923544 0\n",
      "L1\n",
      "output signal MSE: 5.494529043411899\n",
      "output signal SNR: 33.27363612961901\n",
      "iteration: (53,)\n",
      "SCAD 3.5360073640529697 10.98252607049636 3.7\n",
      "SCAD\n",
      "output signal MSE: 0.150944251742575\n",
      "output signal SNR: 48.88477507439896\n",
      "iteration: (79,)\n",
      "MCP 2.773982124634157 2.424422936654266 1.4\n",
      "MCP\n",
      "output signal MSE: 0.6793889731804016\n",
      "output signal SNR: 42.351755925393526\n",
      "iteration: (18,)\n",
      "L1 0.8393417668433006 7.477276129923544 0\n",
      "L1\n",
      "output signal MSE: 5.399118733007213\n",
      "output signal SNR: 33.349712080903615\n",
      "iteration: (51,)\n",
      "SCAD 3.398471498078903 13.405268722005113 3.7\n",
      "SCAD\n",
      "output signal MSE: 0.23174362754618033\n",
      "output signal SNR: 47.022862855041026\n",
      "iteration: (90,)\n",
      "MCP 1.281579539485285 1.4089268310415122 1.4\n",
      "MCP\n",
      "output signal MSE: 1.0740306734062655\n",
      "output signal SNR: 40.36277401692796\n",
      "iteration: (501,)\n",
      "L1 0.8393417668433006 7.477276129923544 0\n",
      "L1\n",
      "output signal MSE: 5.285258470107026\n",
      "output signal SNR: 33.44227855430074\n",
      "iteration: (51,)\n",
      "SCAD 6.481715836209282 9.214141892063944 3.7\n",
      "SCAD\n",
      "output signal MSE: 0.26966077002954725\n",
      "output signal SNR: 46.364763159807026\n",
      "iteration: (35,)\n",
      "MCP 1.644030343266579 1.466933144769021 1.4\n",
      "MCP\n",
      "output signal MSE: 1.1674725632625815\n",
      "output signal SNR: 40.00047403299391\n",
      "iteration: (18,)\n",
      "L1 0.8393417668433006 7.477276129923544 0\n",
      "L1\n",
      "output signal MSE: 5.403949766704615\n",
      "output signal SNR: 33.345827830023254\n",
      "iteration: (52,)\n",
      "SCAD 6.238450208509903 13.084694890109969 3.7\n",
      "SCAD\n",
      "output signal MSE: 0.24592816288442995\n",
      "output signal SNR: 46.764858207535525\n",
      "iteration: (46,)\n",
      "MCP 2.7661286140827364 4.389101518743512 1.4\n",
      "MCP\n",
      "output signal MSE: 0.44504006066155366\n",
      "output signal SNR: 44.18894980186333\n",
      "iteration: (28,)\n",
      "L1 0.8393417668433006 7.477276129923544 0\n",
      "L1\n",
      "output signal MSE: 5.243779326559557\n",
      "output signal SNR: 33.47649679284584\n",
      "iteration: (53,)\n",
      "SCAD 5.967419311084065 11.304418107101268 3.7\n",
      "SCAD\n",
      "output signal MSE: 0.20527751701693386\n",
      "output signal SNR: 47.54952700365628\n",
      "iteration: (40,)\n",
      "MCP 6.458465369621139 11.266950827769213 1.4\n",
      "MCP\n",
      "output signal MSE: 0.3531826786760733\n",
      "output signal SNR: 45.19294690220932\n",
      "iteration: (70,)\n"
     ]
    }
   ],
   "source": [
    "scad_opts\n",
    "mcp_opts\n",
    "\n",
    "idx = 6 # Input SNR = 20dB\n",
    "\n",
    "l1_times = []\n",
    "scad_times = []\n",
    "mcp_times = []\n",
    "l1_snrs = []\n",
    "scad_snrs = []\n",
    "mcp_snrs = []\n",
    "\n",
    "ses = np.linalg.norm(Y - np.tile(y_true.reshape(-1, 1), (1, d)), 2, axis=0) ** 2\n",
    "\n",
    "mse = np.mean(ses)\n",
    "mean_snr = 10*np.log10(y_true_norm_sq/mse)\n",
    "    \n",
    "print ('sigma squared:', sigma_sq)\n",
    "print ('input signal MSE:', mse)\n",
    "print ('input signal SNR:', mean_snr)\n",
    "\n",
    "for trial in range(num_trial):\n",
    "    scad_opt = scad_opts[trial]\n",
    "    mcp_opt = mcp_opts[trial]\n",
    "    Y = Ys[trial]\n",
    "    PARAM_GRID = [\n",
    "            {'penalty_f': [\"L1\"], 'gamma_ratio': [L1_vec[idx,0]], 'rho_ratio': [L1_vec[idx,1]], 'penalty_param': [0]},\n",
    "            {'penalty_f': [\"SCAD\"], 'gamma_ratio': [scad_opt['gamma_mult']], 'rho_ratio': [scad_opt['rho_mult']], 'penalty_param': [3.7]},\n",
    "            {'penalty_f': [\"MCP\"], 'gamma_ratio': [mcp_opt['gamma_mult']], 'rho_ratio': [mcp_opt['rho_mult']], 'penalty_param': [1.4]}\n",
    "        ]\n",
    "    for params in ParameterGrid(PARAM_GRID):\n",
    "        penalty_f = params['penalty_f']\n",
    "        gamma = params['gamma_ratio']*sigma_sq\n",
    "        penalty_param = params['penalty_param']\n",
    "\n",
    "        if penalty_f == 'L1':\n",
    "            l1_beta = None\n",
    "            rho = params['rho_ratio']*gamma\n",
    "        else:\n",
    "            rho = max(params['rho_ratio']*gamma, 1/penalty_param)\n",
    "\n",
    "        print (penalty_f, gamma, rho, penalty_param)\n",
    "\n",
    "        t = time()\n",
    "        invX = V.dot(np.diag(1/(1+rho*S))).dot(V.T)\n",
    "        B, obj, err_path = admm(Y, gamma, rho, Dk, penalty_f, penalty_param,max_iter=500, tol_abs=10 ** (-3), \n",
    "                               tol_rel=10 ** (-2), B_init=l1_beta, invX=invX)\n",
    "        t1 = time()-t\n",
    "        \n",
    "        ses = np.linalg.norm(B - np.tile(y_true.reshape(-1, 1), (1, d)), 2, axis=0) ** 2\n",
    "        mse = np.mean(ses)\n",
    "        snr = 10*np.log10(y_true_norm_sq/mse)\n",
    "\n",
    "        if penalty_f == 'L1':\n",
    "            l1_times.append(t1)\n",
    "            l1_beta = B\n",
    "            l1_snrs.append(snr)\n",
    "        elif penalty_f == 'SCAD':\n",
    "            scad_times.append(t1)\n",
    "            scad_snrs.append(snr)\n",
    "        else:\n",
    "            mcp_times.append(t1)\n",
    "            mcp_snrs.append(snr)\n",
    "\n",
    "        print (penalty_f)\n",
    "        print ('output signal MSE:', mse)\n",
    "        print ('output signal SNR:', snr)\n",
    "        print ('iteration:', np.shape(err_path[3]))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[7.580904245376587,\n",
       " 7.2566752433776855,\n",
       " 6.836250066757202,\n",
       " 6.927603006362915,\n",
       " 6.7369489669799805,\n",
       " 6.966160774230957,\n",
       " 6.606335878372192,\n",
       " 6.747577905654907,\n",
       " 6.9270570278167725,\n",
       " 7.115587949752808]"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "[5.341256380081177,\n",
       " 4.859194040298462,\n",
       " 7.01589298248291,\n",
       " 10.355015993118286,\n",
       " 14.445314645767212,\n",
       " 9.923240184783936,\n",
       " 11.524209976196289,\n",
       " 4.892420291900635,\n",
       " 6.257377862930298,\n",
       " 5.464836835861206]"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "[3.3509981632232666,\n",
       " 4.587734937667847,\n",
       " 3.6941120624542236,\n",
       " 8.091553211212158,\n",
       " 58.67125988006592,\n",
       " 2.905949115753174,\n",
       " 58.19444990158081,\n",
       " 2.963503122329712,\n",
       " 4.051483154296875,\n",
       " 9.111088991165161]"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "[33.13767177991045,\n",
       " 33.3553745887724,\n",
       " 33.083594660839914,\n",
       " 33.10341091020146,\n",
       " 33.25117753970848,\n",
       " 33.27363612961901,\n",
       " 33.349712080903615,\n",
       " 33.44227855430074,\n",
       " 33.345827830023254,\n",
       " 33.47649679284584]"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "[46.78086640516967,\n",
       " 45.97396938544265,\n",
       " 48.14785197700323,\n",
       " 47.70647420821572,\n",
       " 47.1029750276221,\n",
       " 48.88477507439896,\n",
       " 47.022862855041026,\n",
       " 46.364763159807026,\n",
       " 46.764858207535525,\n",
       " 47.54952700365628]"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "[43.90870256285072,\n",
       " 43.904481972752045,\n",
       " 38.764425594627184,\n",
       " 47.03811835646353,\n",
       " 44.843583300628104,\n",
       " 42.351755925393526,\n",
       " 40.36277401692796,\n",
       " 40.00047403299391,\n",
       " 44.18894980186333,\n",
       " 45.19294690220932]"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "l1_times\n",
    "scad_times\n",
    "mcp_times\n",
    "l1_snrs\n",
    "scad_snrs\n",
    "mcp_snrs\n",
    "np.savez('minnesota-runtime'+'.npz', l1_times=l1_times, scad_times=scad_times, mcp_times=mcp_times, \n",
    "         l1_snrs=l1_snrs, scad_snrs=scad_snrs, mcp_snrs=mcp_snrs, scad_opts=scad_opts, mcp_opts=mcp_opts, \n",
    "        Ys=Ys, input_snr=INPUT_SNR)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.collections.PathCollection at 0x11ece1c88>"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.collections.PathCollection at 0x1233c3c18>"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.collections.PathCollection at 0x11ed023c8>"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "Text(0.5, 0, 'Gain in SNR (dB)')"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "Text(0, 0.5, 'Runtime (sec)')"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.legend.Legend at 0x11ece1f60>"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAU0AAAFGCAYAAADuN80HAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4xLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvDW2N/gAAIABJREFUeJzt3XtsXNd9J/Dvb0hKomzZfJh+1I/KwzioZKsbDSkHTdzGjYbdwoiBIB2KTbDYYKOIVLeBXSyyYlwgroEAq5BZoEjcwBnKzq6FFAFFxgvDjtOE4zq2vFZqkSNp1VpGZI7tVK0TyyRHkWRJfMxv/7j3ju4M7zzu8M6d1/cDEOTcmbnnDDn8zjn33HuOqCqIiKg4gUpXgIioljA0iYhcYGgSEbnA0CQicoGhSUTkAkOTiMgFhiYRkQsMTSIiFxiaREQuNFe6Am7dcMMNunnz5kpXg4jqzMzMzAeq2lXocTUXmps3b8b09HSlq0FEdUZE3i3mceyeExG5wNAkInKBoUlE5AJDk4jIBYYmEZELDE0iIhcYmkRELjA0ierMucvncPd378a5y+cqXZW6xNAkqjPP//J5vPHBG/jx6R9Xuip1iaFJNYctqfyePvG08f340xWuSX2qucsoiewtqS9s+0Klq1Nxz5x6Bj9/5+fp26+8+woA4OV3X8ZDP3kovf3+zffjc1s+53f16g5Dk2qOvSXF0ASWVpbwxPQTWE4tZ2y/snIFj7/+OACgOdCM+26/rxLVqzsMTap6bEnlN3DPALbdtA0P/vBBvHf+PVxavpS+r7W5FbdsugXPff45bO3aWsFa1g+GJlU9tqQK29q1FTODM7hh9IaM7Ysri4gPxnH9husrVLP6w4EgqnoD9wzgxN4TCLYH0drcmnFfa3Mrgu1BnNh7Arvu2VWhGlaHw+8exsaWjWgONKNJmtAcaMbGlo04/KvDla5aXWFoUk2wWlKLK4sZ262WFLuewMETB3Fh8QK237wdr+1+Ddtv3o4Lixdw8MTBSletrjA0qWawJZXf6fnTePRTj+LI7iO499Z7cWT3ETz6qUdxeu50patWVxiaVDPYksrv+N7jeOz+x9AUaAIANAWa8Nj9j+HY3mMVrll9YWhSzWBLiqqBqGql6+BKb2+vco0gIvKaiMyoam+hx7GlSUTkAkOTiMgFhiYRkQsMTSIiFxiaREQuMDSJiFxgaBIRucDQJCJygaFJROQCQ5OIyAWGJhGRC77P3C4iIQBBAFDVSXNbBEASQEhVR/2uExFRsSrR0nzEDMugiITMEIWqxgAkrdtERNXI19A0W5RHAUBVR1U1DmAARisTABIAwn7WiYjIDb9bmjsAdJotzH3mtjYA87bHdGY/SUQGRWRaRKbPnj3rRz2JiBxVons+Z7YwrZZnQao6pqq9qtrb1dVV3toREeXhd2jOweiCA0aXfIf5vcPc1mY+hoioKvkdmpMwR85hBORRAOO2bUEAMZ/rRERUNF9DU1UTMEbIIwA6VXXS1lUPA0hat4mIqpHv52mq6pj546TDNiKiqsYrgoiIXGBoEhG5wNAkInKBoUlE5AJDk4jIBYYmEZELDE0iIhcYmkRELjA0iYhcYGgSEbnA0CQicoGhSUTkAkOTiMgFhiYRkQsMTSIiFxiaREQuMDSJiFxgaBIRucDQJCJygaFJROQCQ5OIyAWGJhGRCwxNIiIXGJpERC4wNImIXGBoEhG5wNAkInKBoUlE5AJDk4jIBYYmEZELDE0iIhcYmkRELjA0iYhc8D00RWTE/D5o2xYRkbCI7PO7PkREblSipTkoIrMAEgAgIiEAUNUYgKR1m4ioGlUiNPeoarcZkgAwACBp/pwAEK5AnYiIilKJ0AxmdcXbAMzb7u/MfoKIDIrItIhMnz171pdKEhE58T00VXXUbGV2ikhRrUpVHVPVXlXt7erqKnMNiYhy8zU0zRZjxLw5ByAIo2veYW5rM7cTEVUlv1ua0wCsY5nd5u1xGOEJ83vM4XlERFWh2c/CVDVutjbnAcyqahwARKTX7KonrW1ERNXI19AEjOOTxWwjIqpGvCKIiMgFhiYRkQsMTSIiFxiaREQuMDSJiFxgaBIRucDQJCJygaFJROQCQ5OIyAVXVwSJyKcB9AEI4eokG/Mw5sGcAXBIVX/raQ2JiKpIUS1NEfmyiLwFYBLGRBvHABwCMGb+LAD2AlgQkXER+d0y1ZeIqKIKtjRF5BCA6wH0q+qxAo+9HsZM7DERGVTVl7ypJhFRdcjb0hSRPQDGVfU/FgpMAFDVc+aEwXcB+AuvKklEVC3ytjRV9UCpO1bVXaU+l4ioWrkePReRzSKyOWvbp7O3ERHVo1JOOYri6kzrlnYAI2uvDhFRdSslNHtV9R/tG1T1R+DSu0TUAEoJTRGRTU7b11oZIqJqV0poTiCrKy4iT8A4b5OIqK65XiNIVYdEZEZE5mBcCRQ0v+/0unJERNWmpIXVVLVHRHbCDExVfdHbahERVaeSQtM8vehO6zxOEfkYjPDkdedEVNdch6aIfBnGdebXA3jS3NwN4BEYl1ASVcTS0hLOnDmDy5cvV7oqVKU2bNiA2267DS0tLSXvo5SW5jCAXhizGgEwTjkSEa5dThV15swZbNq0CZs3b4YIT+agTKqKubk5nDlzBnfeeWfJ+yll9LxDVc85bOe7lCrq8uXL6OzsZGCSIxFBZ2fnmnsipYTmiyLyOQBqq8w4gNiaakLkAQYm5ePF+6OU7vkeAC8C6BaRnwLYAeOUo0+vuTZERFWulPM0zwHotZ1yNMpTjoioUZQyy9HHROQ6MygPAdguIl8Vkeu8rx4ReSkej6Ovrw/d3d3o7u7G8PAwAGB0dDS9rb29He3t7enbo6OjAJCxzfpKJpNFl93e3p7z8bFYLF1OtSule34AQD+A38IIzU4Y6wTtAE85ohp1fnkZ4++/j9OXLuGu1lYM3HgjNjWXdBqz58bGxjA4OLjm/cTjcfT392NiYgKhUAgAMDk5CQDYt28f9u3bBwDpIB0ZWT1x2czMDNra2vKWE4vFEA4XP39Pf38/EomEq+dUUikDQd2q+o65tEUYwKdV9U/AWY6oRr2aTOLWI0fwV2+9hdF//Vf81Vtv4dYjR/Cqi1ZUOUWjUU/2s2fPHoyMjKQDEwAikYgn+7aLx+OuHj8xMYGhoSHP61EupYTmvPk9DOBt21VAHLakmnN+eRkPnDyJ8ysruJhKAQAuplI4v7KCB06exIXl5TWX0d/fj1js6skl3d3dAK52ifv6+tLd1rGxsYyu89DQEOLxOHp6etItwLGxMfT09KCnpye932QymX5Mf3//qjokk0nE4/GyhGSjKaX/ERORozAGgfYDgDkoNO1mJyKyT1VHzZ8jAJIAQtY2Ij+Mv/8+UqqO96VUMX72LHbfcsuayhgYGMDExATC4TDi8ThCoRDi8TiOHj2K2dlZJBKJdEBGo1HMzs4CABKJBILBIGKxGGZmjGtJ4vE4JiYm0rd7enrw4osvpu8bGBhw7FbPz88jGMyeO9y9nTuvzssTDocdy6p3pYye7xWRPwOQzBo1L/q3JyJhGOunj4pIyNxvTESCIhJSVXfte6ISnb50Kd3CzHYxlcJbH3645jIikUi6lTg+Po6BgQGMj48jkUigr68v/TjrPotTyI2Pj2d0ZQcGBnDo0CHs2rULbW1t6eOS2To6OpBIJNK3JycnsX//fiQSCSwsLBT9Wl588cVVxzQTiUT62CgATE1NZdwfiUQ8CexqUeosRz/Kur2WU44GAFi/5QSMbj9Dk3xxV2srrgkEHIPzmkAAH9m40ZNyrNZlLBbDyMgIjh49ikceeSSju2wFaz7Zo89zc3PpEOvo6EhvtwZXAOOYYTAYRCgUSg/SRCIRRCIRtLe3r/m1BYPBVWGdK7zrQcElfM2rf1wRkevME9+d7gupqv3qoTZcPU4KGKPxRL4YuPFGBHJcJRIQwUBXlzflDAwgGo2mg826bbG61uPj4xnbAKCtrQ3JZBLJZBL9/f3p5yWTSUxOTjqOOltd+JmZmXQr78CBAxgaGspocZJ7BZfwFZHvicgQgG+q6kv5Hm+eq/nXAP4MQK4lfDtybM+330EAgwBwxx13uH06UU6bmpvxwrZteODkSaRUcTGVwjWBAAIieGHbNlzr0WlHkUgkfboPYLQ8+/v70dPTAwAYGhrC4OAghoaG0gNFQ0NDCIVCCIfD6OnpQTgcRjQaRTweR3d3N9ra2jAyMoJgMFjU+ZKhUCg9Um0F565d/q60bb1ewGihTk1Nob+/H/F4HPPz84jH46u699VGNMdB8IwHGQM13wRwJ4xrzBMwBm7mYLQMgwBC5vcxGFcJve2wn/TxShGZUtU+ERkBMGUe04wACOYbDOrt7dXpaVdjTtQgTp06hS1btpT03AvLyxg/exZvffghPrJxIwa6ujwLzEbj9jxNv+V6n4jIjKr2Fnp+Ue8KVZ0EMGkO2uyEcSJ7t3l3Ekb3+msAYjlmQLIERSQIo7XZYe5vHMZUc4ARupz4g3x3bXPzmkfJyVDNgekFVx+lZiux5EEaM3yt7nabtU8R6TVH1JMcOSeialaR/oeqjsHoxttvExFVvVKuCCIialgMTSIiFxiaREQuMDSJiFwoKTRFZLO5lK91+2OchJiIGkEpM7d/GcAkjKV8Ld0wJicmqimFLu4o5uIPaiyltDSHYZzgnr5g15zAo77PaKW6c+XYU7jy+rdzBqOq4srr38aVY0/5XLPyqeRyF5OTk+nn9fT0rJqs2Jo3NLtcS19fn+OkJu3t7ejp6UFfXx96enrKP6Gxqrr6AjBnfn8ra/u8232V8tXT06NETt54442iH5tKpfTSL/5Wf/u/PqGXfvG3mkqlXN3vp2g06sl+ZmZmNBgM6szMTHrbxMTEqsft27dP9+3bt2p7W1ubLiwsFCxnampq1bbZ2dmM58/Ozurs7OyqumVvsz9/ZGREQ6FQwXpFo1HHx1lyvU8ATGsRGcR1z6khiQjW3/swWrb0Y+nUREaLU80W5tKpCbRs6cf6ex+u6Hrq9bDchbUGkDWNXTAYzJhjc8+ePYhGoxnb7PWMRqOIRCLo7e0tuJzG4OAgOjo6MmbL91IpobkHxkxG3SLyUxGZg3HN+JfzP42ouuQKTq8Dk8tdAL29vYjFYhgeHnacmi4ej+e9Zt2axd4+NV4+9hmlvMZ1z6mhWcEJAEunJrB0yvhH87KFyeUujDlBZ2ZmMDIygp6eHvT29mJiYgJtbW3p15mL/f5wOFzUMctgMFg9oWlz1Pyy5tGEXl1kjahmWMFpBSYAT7vkXO4iku6OR6NRRKNRDA8PY3h4ON0lzzcxcjQaRSwWS/+urHk37d33bIWCeC1ch6aI7AHwvezNMI5xNnlRKSI/WV1yuyuvf9vT4ORyF5mGhoYyDiPY65YtHo+nW9aAEfhW+OYyNTWV8QHkpVKOaX4TxtyZH4E5LyaAdpQwIztRpWUfw7z2i686Dg6tVaMvdxGLxTA2dnUys5GRkYx6O9VtcnLSscUYDodx6NChnGWNjY0hkUiU7fhtKd1zUdVveV4TIp/lGvSxH+MEvOmqN/pyF9YxTKtla72W7Lr19/enX8vQ0BCOHj2acQgDMD5EOjo6MrroO3fuREdHB+bn59Hb25vRMvVaUctdZDxB5JsA/klV/095qpQfl7ugXNwsd1FolLzaTjuqJVzuYrVxADMisgBjraA0Vd1Rwv6IfCcikHWbcgaivcUp6zYxMF2o5sD0QimhOQbjRPbqXjKOqID123dDVXMGohWcDEyyKyU0u1WVgz5UFwoFIgOTspUyej4tIps8rwkRUQ0opaU5BeAdETkEYNZ+h6r+T09qRURUpUoJzQEAb8NY+9w+8KMAGJpEVNdKufa84JA8EVG94hpBREQuFGxpish+AFFVfce8nXMKOFV90ruqERFVn2Jamv0wpoCz7M3xVeY55onK69zlc7j7u3fj3OVzla5K2VRyuYvu7m7HuT6Hh4dXndqVa+mL7KUtipnkxHPFTO9eTV9c7oJycbPchZMfnPiB4jHo3/+/v/eoRt6oh+UuVFWDwaDjMhShUEjb2tpW1dNp6YvsOuSqaz6VWO7CEZfwpVr39Imnje/Hn65wTTLVw3IXlnA4nDH3ZiKRQG9v5thyoaUv7AYGBgouf+G1UpbwPe2wbTuA8kyTTFQmz5x6Bg/95KH01yvvvgIAePndlzO2P3PqmTWVw+Uursqe/i4aja6qb6GlL+z279+/ahakcivlPM3O7A2qekxEeCoS1ZSllSU8Mf0EllPLGduvrFzB468/DgBoDjTjvtvvW1M5XO7iquxZ2q1JmS3FzLhur8PAwEDBCY+9VnRoisjPYJzAfr2I/DTr7iCAeS8rRlRuA/cMYNtN2/DgDx/Ee+ffw6XlS+n7WptbccumW/Dc55/D1q6tayqHy11EMl6L1UUPhUKrXmOhpS9y1cFPblqaEzCWtegDMJl13zy4hC/VoK1dWzEzOIMbRm/I2L64soj4YBzXb7jek3K43MVVQ0ND6VUpnRZJy7f0RTUo+pimqh5Q1TEAk+bP9q8fqbFKJVHNOfzuYWxs2YjmQDOapAnNgWZsbNmIw7867FkZjb7chZ3VmpyamnKse66lL6qF64EgVV3T/PgiEja/RmzbIuY2fw9OEAE4eOIgLixewPabt+O13a9h+83bcWHxAg6eOOhZGZFIBGNjY+mWlX25i56eHkxPTyMUCqWXu+ju7k4P8ljLXQwPDyMcDqfPs9y5c2d6uYti2Je7sMrwa7mLbOFwOGcX2770hVXPSge9XSnLXWwGMAKgDVmLqWmBmdtFJAygX1WHRGQKgNUfCarqpIgMwjhXKuc5BFzugnJxs9yF3ce+9zF89vc+i6//0dfRFGjCSmoF33jlG3j2zWdxbO+xMtS0vlVz1xpY+3IXpYSmlVjj2fepiwXXRGRWVbvNFueUqsbMUA2p6miu5zE0KZdSQ5MaSyXWCArqGmduN7vh1hHgNmSOvK86pYmIqFpUZOZ2syU5JCJFnTcgIoMiMi0i02fPnl1L0UREa+LrzO0iEjIfF4exkuUggCSuHhttAzCX/Txz1H4MMLrnJdSZiMgTfs/cHgZgDfK0ATgK4/xO6zhCEDzfk9ZA86wuSeR2DMeJ3zO3jwHYZY6SQ1UnAUBEes1BoGS+kXOifDZs2IC5uTl0dnYyOGkVVcXc3Bw2bNiwpv2U0tIsmaomYXazs7av2kbk1m233YYzZ86Ax70plw0bNuC2225b0z5ch6Z5ypFjG7fQeZpE5dTS0oI777yz0tWgOldKS3PV+ZkwjnMeXWNdiIiqXinHNFedwC4iB+AcpkREdcWTmdvNY5Vrn6yPiKjKlXJM86sOmzuRdR06EVE9KuWY5p87bEsAqMx0KUREPvL7PE0ioppW9DFNEbku34qTIvJpb6pERFS9igpNEfkejGvEF0Tkf2Tdt9lcP2jK8clERHWkYGiKyB4Y14b3wLjWfJeI7Dbv++8wjmcugANBRNQAijmmOQggoqrvAICI7AXwhIh8zby/R1U5vTURNYRiQjNoBSYAmDOsdwMYVNUny1YzIqIqVMwxTadFkZMMTCJqRMWEptPkHJwImIgaUjHd824RyZ6Mo91hG2c5IqK6V0xoOq0M+aLXFSEiqgUFQ1NVv1boMUREjcKTWY6IiBoFQ5OIyAWGJhGRCwxNIiIXGJpERC4wNImIXGBoEhG5wNAkInKBoUlE5AJDk4g8p5p/Tp9C91czhiYReerKsadw5fVv5wxGVcWV17+NK8ee8rlm3mBoEpFnVBW6eB5LpyYcg9MKzKVTE9DF8zXZ4ixl3XMiIkcigvX3PgwAWDo1AQBYf+/DEJGMwGzZ0p/eXmsYmkTkqVzBWQ+BCTA0iagMsoPTCs9aD0yAxzSJqEzswWmp9cAEKhCaIjJofo3YtkVEJCwi+/yuDxGVh3UM0y7fqHqt8DU0RSQMIKaqYwCCZlCGAGNpYABJ6zYR1a7sQZ9rv/gqWrb05xxVryV+tzSDAMLmzwnz9gCApG1b2OF5RFQjco2Sr7/34boITl8HgswWpiUEYBxAD4B52/bO7OeJyCCAQQC44447yllFIlqDfKcV5TsdqZZUZPTc7ILHVTVezC/MDNsxAOjt7a3NjyeiBiAikHWbco6S24NT1m2qucAEKnfKUVhVh82fkwA6zJ/bAMxVpkpEtFaqivXbd0NVHQPR2l6LLUxLRUbPVXXU/DkMo4seNO8OAoj5XSciWjv7Nee5AtO65rxWAxOozOj5iIjMisgCAKhq3HZf0rpNRLWjEa45t/g9EBQD0O6wfczh4URUIxrhmnMLL6MkIk+4veY8VzfeUuj+SmFoEpFnir3m/Mqxp6CL53O2Oq3WqazbhPXbdxcs188A5rXnROSpQtece3380+9JjxmaROSpQtec57s6yO3xz0oMQLF7TkSecQo96zZwtcXp1ZyblRiAYmgSNYBzl8/hE099Aq/tfg3Xb7i+LGXku+YcWB1qXs256fekxwxNogbw/C+fxxsfvIEfn/4xvrDtC57vv9Rrzq37rO32+9zwc9JjHtMkagBPn3ja+H786bLsv9hrzlu29Gdcc+7lnJt+TXrMliZRHXrm1DP4+Ts/T99+5d1XAAAvv/syHvrJQ+nt92++H5/b8jlPysx3zTlwNdSyA7PQ8c9i5Qpgr4OToUlUh5ZWlvDE9BNYTi1nbL+ycgWPv/44AKA50Iz7br/P03KLGbgB3B//LMTrAM6H3XOiOjRwzwBO7D2BYHsQrc2tGfe1Nrci2B7Eib0nsOueXWWtx7nL53D3d+/Gucvn0tuKOf7pZrJivyc9ZmgS1amtXVsxMziDxZXFjO2LK4uID8axtWtr2etgH4CylHr804nXAVwMds+J6tjhdw9jY8tGXFq+lD7e2NrcisO/OozPfPQzZS/fPgBlH7V3e/wzl0pMeszQJKpjB08cxIXFC+j9nV783QN/h6+88BVM//s0Dp44WJbQ9HIAqtiA8yqAi8XQJKpjp+dP49FPPYqv/9HX0RRowpHdR/CNV76BZ998tizleTUA5XYCjmIHoLwgtTYZaG9vr05PT1e6GkSUwxtn38CDP3wQ751/D5eWL6W3tza34pZNt+C5zz+X93iq1zMgFUtEZlS1t9DjOBBERJ5aywBULcwAz9AkopLlCi1rAKo50IwmaUJzoBkbWzbi8K8O592flzMglQtDk4hKYs1jmUqlVt1nDUD9h2tvwUv3DWP7zdtxYfECDp44WHC/uYKzGgIT4EAQEZXA3o1e+fVxtD74fQQCV9tgv5z/Jf66+0/w1aXfouX8v+Gl0JfxrbsewLNvPlv0II9fE3C4xZYmEbkmIli34yEE2u9CauE0Lj33pXSLU1XxT9v/C4aXz6Ol46NILZxG0/KH+JtP/Q1+8fGvFD3Lul8TcLjF0CSitGIuWbQEAgGjhZkVnFY32tputQ4BuBrksfZl5+XlkKVi95yIAJR2qo8VnJee+xJSC6dx8eAfGtuzAtPaX7GzrK/b8RAWj37Hlwk43GJLk4gKnupjb0Fmn+oTCASw4TOZi5Y5BSZQ3CCPU2CWcwIOt9jSJKoyfixNkS3ftGyXjz2JlV+9mjMIU6kULj+/+iTzdTsecmwN5hvkyRWYheroJ7Y0iaqM08xAfnBqzaVSqXRgBtrvWhWEqVQq3TXH+raM/dkHh3KVZbf+3ocRCAQ8mwGpXBiaRFWm3EtT5JMdnBcP/mE6MFMLp7F49DvpbvGqwLySRMuWflzznw87jqrb5VvmYv323XlbkFYdvbyE0g12z4kqLNfMQFOJKex9fi/WNa3DleUreOGtF/DPf/HPZe+yOy121vrg99PdZsDoejsFphV29sGhS899KeM8Ti9mWa/kaUcMTaIKyzUzkEIRnYkCAAISQEpTZVtNMqNch1bg4tHvYN0OY2o3+3FIp8AEVo+qLx79Tro77uUyF5XA7jl5xmlpAyos39IUAQkg2B7Ex2/9OIDyd9mzW4HXfvHVdFfdHpyW5o9+NufxRys4reOPgHNgAuWbZb0c2NIkz5R7be169uYHbyJ8ZxgH4gcytqc0hdbmVrz+b68DAF5656WyrSZZzGJnK78+nvEcWb6Yc5QcMIIzIxx9nmW9HBia5JlcSxtQYUsrS3jy2JNQrG5d/cvZf7n6uNRSWVaTLLTWzrodD2Hl18fTg0LZxzgLDdxY/J5lvRwqEpoiElLVuO12BEASQEhVRytRJ3Kv3GtrV+J8xWKUo14D9wzg+8e+j58lfgaBOIanpdjJfN3It9aOqmLx6HfSgdl0+33pFiTg/jhkNQ/yFMP30BSRMIAogG7zdggAVDUmIsHsQKXqVe61tau1u1+uev3m4m8w/MlhfOv/fitvaFqT+Xr9QeLUCnS6WscaBa+1ARyv+D4QpKoxAAnbpgEYrUyY28N+14lKU+61tSt5vmI+5arX8b3H8cnbP4lr1l2Tnrw3IMa/aJM0uZrMt1ROxxntLVD79G/W/ZU+2dxv1XBMsw3AvO12Z6UqQu5ZSxvcMHpDxvZSWkPl7u6Xys96Za8e+ac/+FMsXF7Adeuvwz/8p38o+2qSTurhOKSXqiE0CxKRQQCDAHDHHXdUuDaUzau1tcvd3S+Vn/XKXj3y9utux+/f9PtIXkri3lvvLftqkrnU+nFIL1VkNUoRmVLVPvPnEQBT5jHNCIBgvsEgrkZZffoP9eNHp360am3tyNYIDvUfcrWvta5kWC7VWi/yTi2tRjkOIGj+HAQQq2BdqARW6+jI7iPp1tCjn3oUp+dOu97XWlYyLKdqrRf5z/fQNFuTveZ3WCPl5qh6kiPntef43uN47P7H0BRoAgA0BZrw2P2P4djeYyXtr9SVDMutWutF/qrE6Pmkqrar6qRt25iqxlR1zO/6UPWxBkO237wdr+1+zdVKho1YL/JXNXTPiTJ42d1vhHqRvyoyELQWbgeCzi8vY/z993H60iXc1dqKgRtvxKbmmjhpgIh8VOxAUF2nx6vJJB44eRIpVVxMpXBNIID/NjuLF7Ztw31tbYV3QESUpW5D8/zyMh44eRLnV1bS2y6aM0j/8fEYHN7WAAAILUlEQVTj6O/qwk3r1uHua65h65OIila3STH+/vtYtgWm3TKAH549CwDYIMLWJxEVrW4Hgsbffx+XCj8Ml1VxfmUFD5w8iQvLy4WfQEQNrS5D898vX0YsmSz8QJulVArjZuuTiCiXugzNR95+2/VzLqvijYsXy1AbIqondRmab374YUnPm1ta8rgmRFRv6jI0f2/jxpKe18kRdCIqoC5Dc/+dd7p+zgYRbL322jLUhojqSV2G5u9s2IDvfuQjrp7TEghgoKurTDUionpRl6EJAP/1ttvw3h/8Ab54003Y2tqa84VuEMGmpia8sG0brmX3nIgKqOuUuHn9evzvLVsAABeWl/H0b36DH3/wAZYBdLW04OZ167D1mmsw0NXFwCSiojRMUlzb3Iy/vPVW/OWtt1a6KkRUw+q2e05EVA4MTSIiFxiaREQuMDSJiFxgaBIRucDQJCJygaFJRORCzS2sJiJnAbzrY5E3APjAx/JYPstn+ZUp/3dVteC11DUXmn4TkeliVqhj+Syf5ddf+U7YPScicoGhSUTkAkOzsDGWz/JZfsOWvwqPaRIRucCWJhGRCwxNByISyrF9XyXKF5GQiEREJFKh8iMiEhaRQT/Kb2QiMmh+jTjc58v7j/JjaGYRkTCAiRzb+ypU/iOqOgkgmCvQy1W+WV5CVWMAEj6Uvyo0bKFd9tDIUX7OIPO47DCAmKqOwfhbh7Pu8+P95/T6ffvQLvD3r4oPbYZmFiscqqV884161LxvVFXjfpZvst7AwXKW7xQaVkib9UqWM7RzlJ8zyMogCMDaf8K87Zs8r9WXD+08f3/fPrSLwdAsgoiEzD9aJewA0Gl+2vvePTNDMiEiCwDmy1ycU2gMAEjatvkdWr4FmaqOmYEBACEA04Cv779Vr9XnD+1cv2tfPrSLxdAsTkeFy5+z3ix+Hde0iEgbjNDaD+CAiPgdGm3IDOtOP8vPFWTlZLam4raA8OX9l+O1+vahneP37+eHdlEYmgVUuJUJAHO42l1OwngT+2kQwH5VHQWwB4Afx7WyQ8NXTuX7XKewqg5b5fr9/nN4rb5+aNvL9/NDu1gNs7DaGgTNP1QHgA7zTeznP/MkrgZVG8yuUiWo6qRPB+PToQHjH8ZqabXB+BDxs/x82zwnIoPmB5R1jK+tAu8/+2t1+tCe9LF860M7KSIJGP8Lo2UuPy+2NLOYn6S91ieqqk6aB8EB45/W7/ITMAZAIgA6bXXxq/xRAIPmCOagrftUrvKzQ2McV49tBQGUtdXlUL7jtjKVHQYwIiKzZne0Eu+/7Nc6iau//7J/aOf7XZu/h6TjE33EK4KoathOd5qH0bLqV9WY2bpNwBgIKFtoO5Vv3rWqTuWqQyUV+P3PA9hRztZ2nvL3wfj7d5T7Q7sYDE0iIhfYPScicoGhSUTkAkOTiMgFhiYRkQsMTaIyKMf5rNUyYUWjY2g2MBEJisiEeV6gmt9HzKswCj131otZfzzcz76s1xG1vw7zdamIOM1gNWO/RND22FnrnEnzMUWFlohEkeeSP7OuU0WUl331zbxT/clfDM0GZf5DzsI4WblHVQXGeYnFziIzBCDqQVXWvB8zSAYA9Nleh9P12gkAkSIvxYurarf51W7uc8gMxHx1icA4n9TtRQjZ5Q3DuGwwHebmPtvY4qwsnqfZgMwW2AKAoWo4WXitREQBdJtXT+V6zAhsMxSpar/tvhkA47YrUUZgXMrXk7WPEIAZAO2q6nhliojMwvi95jwB3gzCPlXtK1BeBMABM0TtdZhQ1e5c+6fyYkuzMY3AmKOw5gPTptgW8n4U39p0xQy0oIdXDK1qLZvXnXeU83JOyo+h2Zh6UcSkC2LMoG0dY5uVzJnEM47xWccFRWTKdkyuYDB5tJ9RABPmccxIvmOyZujEYHR/i2a7xG80VysTxlyQqwLTPHY8Y76eKRSYqUpErC74CK5eymkXgw+zuJMzhmZjCqK42YLmYRzvbIfxD1xoEOIRGNcLW93JUgd4XO3HvB66H8brmgCwIPnnfhyGMQlJvjAOmSG3YHb/owBGClx7vQPOs/5PwZgbst3skjuFero8GIdORpD7OvcEim9Zk8cYmo0pgazJfK1wML/sMzwlzZ/HYAxC5PtnPWRrhY2j9Fl5XO/HrKt9IGgk14BJka3NuBly7TAGq1DE4Qxr7sc0s4UaVNUh2+YprJYuL+s1OA08zaHyE2M3LIZmY5pG1mTCtn/W9D+92U0cMbuWM0Xsdzbrdqn/2GvajznKPArnrq3Fam0WE8hjQNGrQWa34IMAXM9/qaoxc2Bo0OHUo4pPj9bIGJqNaRjG5Mo5T10xw+RtAEdVtSd7ZLcGrGr12dlam8UeQrBafvm69EmsXo5jresK2SdhtnSjSpZ+aEQMzQZkdn37AUTNQZc2YNV65x0wgqciaxMVS4y1axbMQaugbRBlEMZIeT7D5uOKaW3GYQye5TtPcx5ZAWkdk7S62ebvemj1UzOZr8OaW/JQ1t1tqOCKqY2OodmgzC5sN4zBixlzsGMCxj9ozDzncRLArHnuod9rExXFDLN+82sG5vmnMM6DzNsttrU2i20J7gFgLevrZArOAzQ9MGbDXwBwAM6zz4fsVwTBeC3WQFx2i7nXvJ8qgCe3E3nEdtFAzpPfa6EMyo8tTSKPmCE2BuOUqXLZBWCSgVk5bGkSecgcKJqCc7fai/3Pwjj0wGOaFcKWJpGHzDDrR+ELAVwzB4aGGJiVxZYmEZELbGkSEbnA0CQicoGhSUTkAkOTiMgFhiYRkQsMTSIiF/4/BBRE0uex0rcAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 360x360 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig = plt.figure(figsize=(5,5))\n",
    "plt.scatter(np.array(l1_snrs)-INPUT_SNR, l1_times, label='vector-GTF+L1', s=50, c='c')\n",
    "plt.scatter(np.array(scad_snrs)-INPUT_SNR, np.array(l1_times)+np.array(scad_times), label='vector-GTF+SCAD', marker='x', s=100, c='#f29f3a')\n",
    "plt.scatter(np.array(mcp_snrs)-INPUT_SNR, np.array(l1_times)+np.array(mcp_times), label='vector-GTF+MCP', marker='*',s=100,c='g')\n",
    "plt.xlabel('Gain in SNR (dB)',fontsize=16)\n",
    "plt.ylabel('Runtime (sec)',fontsize=16)\n",
    "plt.legend(fontsize='large')\n",
    "fig.savefig('minnesota-runtime.pdf',format = 'pdf',bbox_inches='tight')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
